Identification of kaempferol as an OSX upregulator by network pharmacology-based analysis of qianggu Capsule for osteoporosis.

Osteoporosis is the most common metabolic disease of skeleton with reduced bone density and weaker bone. Qianggu Capsule as a traditional chinese medicine has been widely used to treat osteoporosis. The potential pharmacological mechanism of its active ingredient Gusuibu is not well understood. The purpose of this work is to analyze the anti-osteoporosis function of Gusuibu based on network pharmacology, and further explore the potential mechanism of Qianggu Capsule. The active compounds and their corresponding targets of Gusuibu were obtained from TCMSP, TCMID, and BATMAN-TCM databases. Potential therapeutic targets for osteoporosis were obtained through DisGeNET, TTD, GeneCards, MalaCards, CTD, and OMIM databases. The overlapping targets of Gusuibu and osteoporosis were obtained. GO and KEGG pathway enrichment analysis were performed. The "Gusuibu-active compounds-target genes-osteoporosis" network and protein-protein interaction (PPI) network were constructed, and the top hub genes were screened by using the plug-in CytoHubba. Molecular docking was used to verify the binding activity of hub genes and key compounds. We identified 21 active compounds and 140 potential therapeutic targets that may be related to Gusuibu and 10 hub genes (AKT1, IL6, JUN, TNF, MAPK3, VEGFA, EGFR, MAPK1, CASP3, PTGS2). Molecular docking analysis demonstrated that four key active small molecules in Gusuibu (including Luteolin, Naringenin, Kaempferol, and Beta-sitosterol) have excellent binding affinity to the target proteins encoded by the top 10 hub genes. Our new findings indicated that one key active compound kaempferol activated the expression of osteoblast specific transcription factor OSX through JNK kinase pathway.

Introduction

Osteoporosis is the most common metabolic disease of skeleton with reduced bone density and weaker bone structure. Patients are more likely to have fractures with the pain and other complications, and the quality of life is dramatically affected (Zheng et al., 2020). Postmenopausal osteoporosis is one common form of primary osteoporosis. Osteoporosis is classified into the categories of “Gu Bi” and “Gu Lou” in the “Huang Di Nei Jing” of Traditional Chinese Medicine (TCM) theory (CSoOaBM, 2019). Bone metabolism is a continuous process, and bone mass will be gradually lost after the age of 35. People especially postmenopausal women and the elderly should pay more attention to osteoporosis prevention (Li et al., 2020a; Saul and Kosinsky, 2021). The fragility fractures induced by osteoporosis will lead to increased disability and mortality, resulting in a heavy family, social and economic burden (Kanis, 2002; Xiong et al., 2021). The commonly used anti-osteoporosis drugs include estrogen, calcium, vitamin D, bisphosphonates, and denosumab; however, the consequent adverse effects are still a challenge to be solved due to the increased duration of use and dose (Jackson and Mysiw, 2014). Therefore, it is needed to explore potential anti-osteoporosis drugs with high efficacy and safety.

It is well known that TCM has a long history of practice in China and other Asian countries to treat a wide range of diseases, including osteoporosis (Li et al., 2017a; Wei et al., 2017; Wu et al., 2017; Zhang et al., 2017; Ge et al., 2018; Li et al., 2020a; Li et al., 2020b; Liu et al., 2020; Zheng et al., 2020). Among those to treat osteoporosis, the clinically commonly used Traditional Chinese Medicine prescriptions include Qianggu Capsule, XianlingGubao Capsule, DuhuoJisheng Decoction, LiuweiDihuang Pill, and Erxian Decoction (Li et al., 2017a; Wei et al., 2017; Wu et al., 2017; Ge et al., 2018; Li et al., 2020a). Qianggu Capsule (Drug approval number: Z20030007, Qi-Huang Pharmaceutical CO. LTD., Beijing, China) is a China Food and Drug Administration (CFDA)-approved TCM for the treatment of osteoporosis and bone loss. The active ingredients in Qianggu Capsule are total flavonoids extracted from Rhizoma Drynariae (English name is Fortune’s Drynaria Rhizome, and Chinese name is Gusuibu), the dried rhizome of Drynaria fortune (Kunze) J. SM. (Wang et al., 2008; CSoOaBM, 2019). The active ingredient of Qianggu Capsules is the herbal Gusuibu. Gusuibu is still a kind of raw material, and there are many compounds in Gusuibu, including flavonoids. Studies have shown that the flavonoids in Gusuibu may treat osteoporosis by improving bone density and reducing bone loss (Zhang et al., 2017; Shen et al., 2020). However, the mechanism of function of many TCM cannot be elucidated due to the complexity of the ingredients in TCM.

In recent years, many researchers have started to use bioinformatics and network pharmacology to analyze the multi-component, multi-target and multi-pathway characteristics of TCM in order to elucidate mechanisms of its action, thus providing directions for further research (Shuai et al., 2020). Gusuibu contains a variety of components. Although some studies have suggested that total flavonoids play a key role in anti-osteoporosis (Wang et al., 2008; Zhang et al., 2017), the molecular mechanism has not been elucidated yet.

In this study we used a network pharmacology approach to analyze the key genes, active compounds and pathways in the anti-osteoporosis of Gusuibu, and further explored possible underlying molecular mechanisms. We identified one key active compound kaempferol in Gusuibu as an upregulator of osteoblast specific transcription factor OSX.

Materials and methods

Exploring potential pharmacodynamic compounds and relevant targets of gusuibu

Using TCM Systems Pharmacology (TCMSP, Version: 2.3, https://tcmspw.com/tcmsp.php) database (Ru et al., 2014), BATMAN-TCM platform (http://bionet.ncpsb.org/batman-tcm/) (Liu et al., 2016) and TCM Integrated Database (TCMID, http://www.megabionet.org/tcmid/) (Huang et al., 2018), the corresponding compounds of Gusuibu and related information were obtained. The compounds retrieved in the previous step were screened for active compounds according to the absorption, distribution, metabolism and excretion (ADME) protocol (oral bioavailability (OB) ≥ 30 and drug-likeness (DL) ≥ 0.18) (Tsaioun et al., 2016; Sheng et al., 2020). The databases were then used to mine potential targets corresponding to the active compounds. The target proteins were then converted to the corresponding gene names and UniProt ID using the UniProt database (https://www.uniprot.org/) for “Homo sapiens” (UniProt Consortium, 2018).

Exploring the relevant targets of osteoporosis

The targets related to osteoporosis were obtained through retrieving GeneCards (https:// www.genecards.org/) (Relevance score ≥10) (Stelzer et al., 2016), MalaCards (https://www.malacards.org/) (Rappaport et al., 2017), DisGeNet database (https://www.disgenet.org/, v7.0) (Score ≥ 0.01) (Piñero et al., 2020), Therapeutic Target Database (TTD) (http://db.idrblab.net/ttd/, Last update by 1 June 2020) (Wang et al., 2020), Comparative Toxicogenomics Database (CTD) (http://ctdbase.org/, Last update by June 2020) (Inference score ≥ 15) (Davis et al., 2020), and Online Mendelian Inheritance in Man (OMIM) (https://omim.org/, updated 25 November 2020) (Amberger and Hamosh, 2017) using the keyword “osteoporosis or postmenopausal osteoporosis”. The potential targets obtained from the these 6 databases were integrated and deduplicated to generate the osteoporosis-related target set.

Network construction and protein-protein interaction analysis

Gusuibu-related targets and osteoporosis-related targets were input to the Vennonline tool (http://www.bioinformatics.com.cn/) to obtain the common targets. These were the candidate targets of Gusuibu to treat osteoporosis. Cytoscape software (version 3.7.2) was used to construct the network diagram of “Gusuibu-active compounds-target genes-osteoporosis” (Shannon et al., 2003).

PPI is the basis of most biological processes and is important to understand cell physiology in both normal and disease states (Xiao et al., 2020). The STRING database (http://string-db.org/; version 11) was used to perform a PPI network analysis on the common targets (Szklarczyk et al., 2017). The species is limited to “Homo sapiens” with a confidence level greater than 0.4. Cytoscape software (version 3.7.2) was used to construct the PPI network, and the plug-in 12 CytoHubba algorithms [Degree, Maximal Clique Centrality (MCC), Clustering Coefficient, Density of Maximum Neighborhood Component (DMNC), BottleNeck, Maximum Neighborhood Component (MNC), Radiality, Edge Percolated Component (EPC), EcCentricity, Closeness, Betweenness, Stress] were used to find the top 10 hub genes (Chin et al., 2014; Xu et al., 2020).

Enrichment analysis of GO and KEGG pathway

The cluster Profiler package in R (R 4.0.2 for Windows) was used to perform Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) analysis to explore biological processes and signaling pathways in Gusuibu for osteoporosis (Kanehisa et al., 2008; Yu et al., 2012; Consortium, 2015). An adjusted p-value of less than 0.05 was used to determine the enrichment term.

Implementation of molecular docking

The Sankey diagram (http://sankeymatic.com/) reveals the correspondence between herbs, components, and targets by constructing the interrelationship between the active compounds of Gusuibu and the top 10 hub genes. AutoDock Vina was used to perform molecular docking between the top 10 hub genes and key active ingredients to predict their binding affinity (Trott and Olson, 2010). The PDB format of the target proteins and the MOL2 format of active compounds could be obtained from the RCSB protein data (http://www.rcsb.org/) and the PubChem database (https://pubchem.ncbi.nlm.nih.gov/) (Kanis, 2002). The smaller the molecular docking score, the more stable the binding affinity between target proteins and active compounds (Chang et al., 2020).

Cell culture

C2C12 mesenchymal stem cells (ATCC) were cultured in DMEM with 10% fetal bovine serum (FBS, Gibco), 100U/mL penicillin G sodium and 100 μg/ml streptomycin sulfate in a humidified atmosphere at 37°C. Cells were plated in 6-well plates, cultured to 60–80% confluence and treated with kaempferol.

Reverse transcription-quantitative PCR

Kaempferol (1,000 ng/ml; Sigma) was used to stimulate C2C12 cells for 24 h before harvest. TRIzol (Invitrogen Life Technologies, Carlsbad, CA, United States) was used to extract total RNA from cultured cells as previously described (Tang et al., 2012). RNA purity and integrity were determined by the RNA 6000 Nano assay with an Agilent Bioanalyzer 2,100 (Agilent Technologies, Santa Clara, CA). 1 μg total RNA was reverse transcribed into cDNA using the GoScript Reverse Transcription System (Promega Corp., Madison, WI, United States) based on the manufacturer protocol. qPCR was performed in triplicate using SYBR-Green SuperReal PreMix Plus [Tiangen Biotech (Beijing) Co., Ltd., China] and the iQ5 PCR system (Bio-Rad Laboratories, Inc., Hercules, CA, United States). The reaction conditions were as follows: 95°C for 30 s, and 40 cycles of 95°C for 10 s and 60°C for 30 s. Data were reported as cycle threshold (Ct) values. The 2−ΔΔCt method was used to compare the RNA expressions. RNA levels were normalized to glyceraldehyde-3-phosphate dehydrogenase (GAPDH) levels.

Statistical analysis

All qPCR experiments were performed in triplicate. Data were expressed as the mean ± SD. Comparisons were made between groups by Student’s t test with p < 0.05 being considered as statistically significant.

Results

Pharmacodynamic compounds and potential target genes of gusuibu

A total of 131 compounds in Gusuibu (Rhizoma Drynariae) were obtained through the TCMSP database (Number: 71), TCMID database (Number: 53), and BATMAN-TCM platform (Number: 7). Based on the respective standards of TCMSP database (OB ≥ 30% and DL ≥ 0.18) and BATMAN-TCM platform (Score cutoff>20, adjusted p-value < 0.05) along with the existence of target proteins, a total of 24 active compounds of Gusuibu were screened. The corresponding target proteins were deduplicated, and a total of 331 target proteins were obtained. The UniProt database was used to convert the target proteins predicted by Gusuibu’s bioactive compounds into gene names. Table 1 includes the basic information of active compounds in Gusuibu.

TABLE 1

Basic information of active compounds in Gusuibu.

Molecule ID	Molecule name	OB (%)	DL	2D structure	InChI key	PubChem CID
MOL000006	Luteolin	36.16	0.25		IQPNAANSBPBGFQ-UHFFFAOYSA-N	5280445
MOL000422	Kaempferol	41.88	0.24		IYRMWMYZSQPJKC-UHFFFAOYSA-N	5280863
MOL004328	Naringenin	59.29	0.21		FTVWIRXFELQLPI-ZDUSSCGKSA-N	439246
MOL000358	Beta-sitosterol	36.91	0.75		KZJWDPNRJALLNS-VJSFXXLFSA-N	222284
MOL000449	Stigmasterol	43.83	0.76		HCXVJBMSMIARIN-PHZDYDNGSA-N	5280794
MOL000492	(+)-catechin	54.83	0.24		PFTAWBLQPZVEMU-DZGCQCFKSA-N	9064
MOL000569	Digallate	61.85	0.26		COVFEVWNJUOYRL-UHFFFAOYSA-N	341
MOL001040	(2R)-5,7-dihydroxy-2-(4-hydroxyphenyl)chroman-4-one	42.36	0.21		FTVWIRXFELQLPI-CYBMUJFWSA-N	667495
MOL001978	Aureusidin	53.42	0.24		WBEFUVAYFSOUEA-PQMHYQBVSA-N	5281220
MOL002914	Eriodyctiol (flavanone)	41.35	0.24		SBHXYTNGIZCORC-CYBMUJFWSA-N	373261
MOL005190	Eriodictyol	71.79	0.24		SBHXYTNGIZCORC-ZDUSSCGKSA-N	440735
MOL009061	22-Stigmasten-3-one	39.25	0.76		RTLUSWHIKFIQFU-ZBWVUXHASA-N	91692436
MOL009075	Cycloartenone	40.57	0.79		NAJCQAAOHKVCES-ZZOQNIIDSA-N	12305360
MOL009078	DavalliosideA_qt	62.65	0.51		ZGIKMQJUZZSYOY-QEIWDELWSA-N	NR
MOL009091	Xanthogalenol	41.08	0.32		ALGFNVZQNNGHPA-YRNVUSSQSA-N	14309735
	Chitranone				ITGPISXKMZIRAV-UHFFFAOYSA-N	633072
	Dihydrocaffeic acid				DZAUWHJDUNRCTF-UHFFFAOYSA-N	348154
	Kurarinone				LTTQKYMNTNISSZ-MWTRTKDXSA-N	11982640
	Catechol				YCIMNLLNPGFGHC-UHFFFAOYSA-N	289
	Naringin				DFPMSGMNTNDNHN-ZPHOTFPESA-N	442428
	Campesterol				SGNBVLSWZMBQTH-PODYLUTMSA-N	173183
	Gamma-Sitosterol				KZJWDPNRJALLNS-FBZNIEFRSA-N	457801
	Narirutin				HXTFHSYLYXVTHC-AJHDJQPGSA-N	442431
	Hesperidin				QUQPHWDTPGMPEX-QJBIFVCTSA-N	10621

Abbreviations: OB, Oral bioavailability; DL, Drug-likeness; NR, not reported.

Potential therapeutic targets of gusuibu in the treatment of osteoporosis

By searching the GeneCard, TTD, MalaCard, CTD, DisGeNET, and OMIM databases, a total of 1,578 potential therapeutic targets for osteoporosis were obtained. By constructing the Venn diagram of the targets regulated by Gusuibu’s active ingredients and the potential targets of osteoporosis, a total of 140 intersection targets were obtained, which are potential therapeutic targets of Gusuibu against osteoporosis as shown in Figure 1. Cytoscape software was used to construct the “Gusuibu-active compounds-target genes-osteoporosis” network for Gusuibu against osteoporosis as shown in Figure 2. There are three compounds (Cycloartenone, Narirutin, Naringin) corresponding to targets that are not intersecting genes, so there are 21 active compounds in the network.

FIGURE 1

Venn diagram for the integrated analysis of the related targets of Gusuibu and osteoporosis.

FIGURE 2

“Gusuibu-active compounds-target genes-osteoporosis” network. The red octagon represents osteoporosis; the blue diamond represents Gusuibu; the cyan rectangle represents the potential target; the brown ellipse represents the active compound contained in Gusuibu. The line between two nodes indicates that there is a relationship, and the size of each node indicates the number of relationships.

PPI network construction and exploration of potential hub genes for gusuibu against osteoporosis

A total of 140 potential target genes of Gusuibu for treating osteoporosis were entered into the STRING database to obtain a PPI network, which involved 140 nodes and 2,199 edges (Figure 3). The data obtained was imported into Cytoscape software (version 3.7.2) for further visualization (Figure 4A). According to the 12 CytoHubba algorithms of Cytoscape software, the top 10 hub genes of Gusuibu for the treatment of osteoporosis were screened. A total of 37 different genes and the number of algorithms to which these genes belong were sorted in Table 2. The top 10 hub genes were consistent with the results of the Degree algorithm (Figure 4B). The Sankey diagram was constructed by using the top 10 ranked hub genes (AKT1, IL6, JUN, TNF, MAPK3, VEGFA, EGFR, MAPK1, CASP3, PTGS2) and the corresponding active compounds of Gusuibu. Luteolin targeted most of the hub genes while PTGS2 targeted most of the active compounds (Figure 5).

FIGURE 3

PPI network based on STRING database.

FIGURE 4

PPI network of potential target genes and top 10 hub genes for Gusuibu against osteoporosis. (A) PPI network constructed by using Cytoscape software. (B) The top 10 hub genes was identified by the Degree algorithm.

TABLE 2

Basic information of top 10 hub genes by 12 CytoHubba algorithms.

UniProt ID	Gene symbol	Protein names	Algorithms
P31749	AKT1	RAC-alpha serine/threonine-protein kinase	A, B, C, D, E, G, H, J, K, L
P05231	IL6	Interleukin-6	A, B, C, D, E, G, J, K, L
P05412	JUN	Transcription factor AP-1	A, B, C, D, G, H, J, K, L
P01375	TNF	Tumor necrosis factor	A, B, C, D, E, G, J, K, L
P42574	CASP3	Caspase-3	A, B, D, E, G, J, K, L
P28482	MAPK1	Mitogen-activated protein kinase 1	A, B, C, D, G, J, K, L
P27361	MAPK3	Mitogen-activated protein kinase 3	A, B, C, D, G, J, K, L
P35354	PTGS2	Prostaglandin G/H synthase 2	A, B, C, D, E, G, J, K
P15692	VEGFA	Vascular endothelial growth factor A	A, B, D, E, G, J, K
P00533	EGFR	Epidermal growth factor receptor	A, C, D, J, K, L
P03372	ESR1	Estrogen receptor	C, E, G, L
P14780	MMP9	Matrix metalloproteinase-9	B, E, H
P10275	AR	Androgen receptor	E, H
P04040	CAT	Catalase	C, L
P29965	CD40LG	CD40 ligand	I, H
P03956	MMP1	Interstitial collagenase	I, F
P16581	SELE	E-selectin	I, F
P98170	XIAP	E3 ubiquitin-protein ligase XIAP	I, H
P09917	ALOX5	Polyunsaturated fatty acid 5-lipoxygenase	H
Q92934	BAD	Bcl2-associated agonist of cell death	F
P55210	CASP7	Caspase-7	H
Q14790	CASP8	Caspase-8	I
P20309	CHRM3	Muscarinic acetylcholine receptor M3	F
P05108	CYP11A1	Cholesterol side-chain cleavage enzyme, mitochondrial	H
Q92819	HAS2	Hyaluronan synthase 2	F
P05362	ICAM1	Intercellular adhesion molecule 1	I
P22301	IL10	Interleukin-10	E
P60568	IL2	Interleukin-2	I
Q07820	MCL1	Induced myeloid leukemia cell differentiation protein Mcl-1	I
P10242	MYB	Transcriptional activator Myb	F
P25963	NFKBIA	NF-kappa-B inhibitor alpha	I
Q9BZD4	NUF2	Kinetochore protein Nuf2	F
P00749	PLAU	Urokinase-type plasminogen activator	H
O14684	PTGES	Prostaglandin E synthase	F
P03973	SLPI	Antileukoproteinase	F
P09486	SPARC	SPARC	F
P01137	TGFB1	Transforming growth factor beta-1 proprotein	I

A: Degree, B: Maximal Clique Centrality (MCC), C: Betweenness, D: closeness, E:BottleNeck, F: ClusteringCoefficient, G: Edge Percolated Component (EPC), H: EcCentricity, I: Density of Maximum Neighborhood Component (DMNC), J: Maximum Neighborhood Component (MNC), K: radiality, L: stress.

FIGURE 5

Sankey diagram. The blue vertical bar on the left represents Gusuibu, the middle bars represent the active ingredients of Gusuibu, and the right bars represent the top 10 hub genes. The size of the bars represents the number of correlations.

The Cluster Profiler package in R was used to perform GO enrichment analysis on 140 potential targets of Gusuibu in the treatment of osteoporosis. 2601 GO items (adjusted, p < 0.05) were obtained, including 2,398 biological process (BP) items, 75 cellular component (CC) terms and 128 molecular function (MF) items. As shown in Figure 6, the top 10 enrichment results of GO-BP, GO-CC, and GO-MF were visualized by bubble chart.

FIGURE 6

Bubble chart of the top 10 GO enriched items.

A total of 166 KEGG enriched pathways (adjusted, p < 0.05) were obtained through the above R package. As shown in Figure 7, the top 20 KEGG pathway enrichment terms were visualized by bubble chart.

FIGURE 7

Bubble chart of the top 20 KEGG pathway enriched items.

We searched for osteoporosis in the KEGG database, and compared the retrieved relevant pathways with the enriched 166 pathways. A total of 50 pathways that may be related to osteoporosis were obtained, and we used Cytoscape software to construct a network diagram of potential target genes and pathways as shown in Figure 8.

FIGURE 8

“Target genes-pathways” network. The cyan rectangle represents the potential target genes of Gusuibu against osteoporosis; the light blue ellipse represents the top 10 hub genes; the brown ellipse represents potential osteoporosis-related pathways.

Molecular docking

Based on the Sankey diagram results, molecular docking was performed between the top 10 target proteins (AKT1, IL6, JUN, TNF, MAPK3, VEGFA, EGFR, MAPK1, CASP3, PTGS2) and four key active compounds (Luteolin, Naringenin, Kaempferol, and Beta-sitosterol) using AutoDock Vina. As shown in Figure 9, a heatmap was used to visualize the strongest affinity docking scores of top 10 target proteins and 4 key active compounds. The binding energy between the top 10 target proteins and the active compound is about −6.4 to −10.8 kcal mol−1, indicating that the key active compounds in Gusuibu bind well to the top 10 target proteins.

FIGURE 9

Heatmap of molecular docking scores.

Kaempferol induced the expression of osteoblast genes

To explore the possible role of one key active compound kaempferol in osteoblast genes, we examined the effect of kaempferol on osteoblast differentiation related gene expressions. C2C12 cells were cultured and treated as indicated with kaempferol (1,000 ng/ml, Sigma) for 24 h. Total RNA was isolated and measured by real time RT-PCR. As shown in Figure 10, after kaempferol stimulation osteoblast specific transcription factor Osterix (OSX) expression was upregulated by about 6.1 folds, while the expression of OSX downstream osteoblastic gene osteocalcin (OCN) increased by 5.9 fold. These observations demonstrated that kaempferol induced OSX gene expression.

FIGURE 10

Kaempferol induced the expression of osteoblast genes. C2C12 cells were treated with 1,000 ng/ml kaempferol for 24 h. Total RNA was isolated and measured by real time RT-PCR. The RNA level from the control group was normalized to a value of 1. Values were presented as the mean ± S.D. *: A star indicates statistical significance compared to control group with p < 0.05. (A) Effect of kaempferol on OSX expression; (B) Effect of kaempferol on OCN expression.

We asked which pathway was involved in kaempferol induced OSX gene expression. To explore molecular mechanisms of kaempferol effect on OSX expression, we used a loss-of-function approach to examine possible pathways involved. The selective inhibitors were used: SB203580 as a specific inhibitor for p38 MAPK pathway, and SP600125 is a specific inhibitor for JNK kinase pathway. C2C12 cells were treated with kaempferol. Different inhibitors were added in the culture medium. As shown in Figure 11, kaempferol treatment led to OSX expression increase by 10.8-fold. Addition of 50 μM SP600125 almost abolished the increment of OSX expression induced by kaempferol. Kaempferol-induced OSX activation remained unchanged after treatment with 10 μM SB203580. The data indicate that kaempferol-induced OSX activation is mediated through JNK kinase pathway.

FIGURE 11

Kaempferol-induced OSX activation is mediated through JNK kinase pathway. Gene expressions were determined via qRT-PCR with or without kaempferol stimulation. The RNA level from the control group was normalized to a value of 1. *: A star indicates statistical significance compared to control group with p < 0.05. Specific inhibitors (50 μM SP600125 or 10 μM SB203580) were added as indicated.

Discussion

Osteoporosis has become an important health issue with aging (Zheng et al., 2020). For example, postmenopausal osteoporosis as a very common form of primary osteoporosis affects most postmenopausal women (Jackson and Mysiw, 2014). Due to severe side effects of current clinical osteoporosis drugs, more and more attention has been paid to possible drugs for the prevention and treatment of osteoporosis from natural products (Li et al., 2020a). Qianggu Capsule is a proprietary Traditional Chinese Medicine approved by CFDA to treat osteoporosis, and its active ingredient is total flavonoids of Rhizoma Drynariae (Gusuibu) (CSoOaBM, 2019). Studies have shown that Qianggu Capsule may prevent and treat osteoporosis by improving bone density and reducing bone loss (Wei et al., 2017). The raw material of Qianggu Capsules is the flavonoids of the herbal Gusuibu, but the mixture is not a single compound at all, and the underlying mechanism is not clear yet. Based on the network pharmacological approaches such as databases mining, PPI network construction, GO and KEGG enrichment analysis and molecular docking validation, this study aims to analyze the active ingredients and potential targets of Gusuibu for osteoporosis treatment and to explore possible mechanisms of Qianggu Capsule active ingredients.

Based on the principle of ADME (setting OB ≥ 30 and DL ≥ 0.18) (Tsaioun et al., 2016), we searched three herbal databases, including TCMSP, TCMID and BATMAN-TCM platform, and screened 24 active compounds from Gusuibu (131 initially, after removing duplicates and null targets), and finally 331 corresponding target proteins were found. The protein name of the target protein obtained was converted into a gene ID through the UniProt database. Subsequently, we searched 6 disease-related databases (GeneCards, MalaCards, DisGeNET, TTD, CTD, and OMIM) and obtained a total of 1,578 potential therapeutic targets for osteoporosis. After constructing Venn diagrams with potential targets of Gusuibu, we obtained 140 potential target genes of Gusuibu for osteoporosis treatment. “Gusuibu-active compounds-target genes-osteoporosis” network was established, which involved 163 nodes and 432 interactions. The PPI network of potential therapeutic targets had also been constructed, involving 140 nodes and 2,199 edges. Based on the CytoHubba 12 algorithms of Cytoscape software (Degree, Clustering Coefficient, DMNC, Bottle Neck, MCC, MNC, Radiality, EPC, EcCentricity, Closeness, Betweenness, Stress), the top 10 hub genes were obtained, including AKT1, IL6, JUN, TNF, MAPK3, VEGFA, EGFR, MAPK1, CASP3, PTGS2. The physiological processes regulated by the proteins encoded by these core target genes mainly include inflammatory response, cell proliferation and differentiation, apoptosis, migration, cell cycle progression, angiogenesis, immune response, endocrine metabolism, bone metabolism, reproductive function, growth and nutrition. Therefore, based on these results, it may be speculated that Gusuibu may play a role in resisting osteoporosis through the multiple processes.

Sankey diagram showed the one-to-one correspondence between the top 10 hub genes and the corresponding active compounds contained in Gusuibu [Luteolin, Naringenin, Kaempferol, Beta-sitosterol, (+)-catechin, (R)-naringenin, Aureusidin, Digallate, Eriodyctiol (flavanone), Eriodictyol, DavalliosideA_qt, Dihydrocaffeic acid, Kurarinone, Chitranone, Stigmasterol, Xanthogalenol]. Among them, Luteolin targeted nine hub genes, Naringenin targeted five hub genes, Kaempferol targeted five hub genes, Beta-sitosterol targeted three hub genes, and the others targeted one hub gene each. Luteolin is a flavonoid found in many herbal extracts with antioxidant, anti-inflammatory and anti-tumor activities, and has been shown to reduce glucocorticoid-induced osteoporosis by modulating the ERK/Lrp-5/GSK-3β signaling pathway in vitro and in vivo (Jing et al., 2019). We summarized the relevant literature exploring the possible mechanisms of Qianggu Capsule and its compounds for osteoporosis treatment (Liu et al., 2001; Wong et al., 2007; Wang et al., 2008; Chen et al., 2009; Hung et al., 2010; Guo et al., 2011; Huang et al., 2011; Ouyang et al., 2012; Li et al., 2017b; Wei et al., 2017; Zhang et al., 2017; Jing et al., 2019; Dong et al., 2020; Ruangsuriya et al., 2020; Shen et al., 2020; Liu et al., 2021) in Table 3.

TABLE 3

Studies related to the regulation of bone metabolism by Qianggu Capsule and its active compounds.

Chinese herbal medicine or compounds	Study	Organism	Possible pharmacological mechanisms	References
Qianggu Capsule (total flavonoids of Rhizoma Drynariae)	Wei et al. (2017)	Human	Improving BMD	Wei et al. (2017)
	Zhang et al. (2017)	Human	Improving BMD	Zhang et al. (2017)
	Li et al. (2017b)	Rat	Enhancing the expression of TGF-β1 and promoting bone metabolism	Li et al. (2017b)
	Ouyang et al. (2012)	Human	Improving BMD	Ouyang et al. (2012)
	Huang et al. (2011)	Rat	Promoting the proliferation and decrease the apoptosis of osteoblasts by improving the ratio of Bcl-2 mRNA to Bax mRNA	Huang et al. (2011)
	Chen et al. (2009)	Rat	Promoting bone metabolism	Chen et al. (2009)
Gusuibu (RhizomaDrynariae)	Shen et al. (2020)	Rat	Enhancing angiogenic-osteogenic coupling during distraction osteogenesis by promoting type H vessel formation through PDGF-BB/PDGFR-β instead of HIF-1α/VEGF axis	Shen et al. (2020)
	Hung et al. (2010)	Rat	Promoting osteoblast maturation by regulating bone differentiation-related gene expression	Hung et al. (2010)
	Wang et al. (2008)	Rat	Promoting bone healing	Wang et al. (2008)
	Wong et al. (2007)	Mouse	Promoting bone healing	Wong et al. (2007)
	Liu et al. (2001)	Rat	Promoting bone healing	Liu et al. (2001)
Luteolin	Jing et al. (2019)	Mouse	Luteolin reduced glucocorticoid-induced osteoporosis by modulating the ERK/Lrp-5/GSK-3β signaling pathway	Jing et al. (2019)
Naringin	Dong et al. (2020)	Rabbit	Naringin may be good natural BMP regulator in bone tissue engineering	Dong et al. (2020)
Naringin and Naringenin	Guo et al. (2011)	Mouse	Naringin and Naringenin revealed a double directional adjusting function of estrogenic and anti-estrogenic activities	Guo et al. (2011)
Kaempferol	Liu et al. (2021)	Rat	Kaempferol promotes BMSC osteogenic differentiation and improves osteoporosis by downregulating miR-10a-3p and upregulating CXCL12	Liu et al. (2021)
Beta-sitosterol	Ruangsuriya et al. (2020)	Human	Beta-sitosterol triggered the molecules involved in bone formation and the molecules inhibiting bone resorption	Ruangsuriya et al. (2020)

Abbreviations: BMD, bone mineral density; BMSC, bone marrow mesenchymal stem cell.

We further analyzed the potential therapeutic genes of Gusuibu against osteoporosis by GO and KEGG pathway enrichment analysis. We used osteoporosis as a search term in the KEGG database and obtained ten directly related pathways, including osteoclast differentiation (hsa04380), AGE-RAGE signaling pathway in diabetic complications (hsa04933), Endocrine and other factor-regulated calcium reabsorption (hsa04961), Mineral absorption (hsa04978), Wnt pathway (hsa04310), MAPK pathway (hsa04010), Apoptosis (hsa04210), Chemokine pathway (hsa04062), T cell receptor pathway (hsa04660), and B cell receptor pathway (hsa04662). Five of the top 10 hub genes (AKT1, TNF, MAPK3, JUN, MAPK1) are enriched in osteoclast differentiation. Eight of the top 10 hub genes (AKT1, IL6, CASP3, TNF, VEGFA, MAPK3, JUN, MAPK1) are enriched in AGE-RAGE pathway in diabetic complications. Eight of the top 10 hub genes (AKT1, CASP3, TNF, VEGFA, MAPK3, JUN, EGFR, MAPK1) are enriched in MAPK pathway. Six of the top 10 hub genes (AKT1, CASP3, TNF, MAPK3, JUN, MAPK1) are enriched in Apoptosis. Three of the top 10 hub genes (AKT1, MAPK3, MAPK1) are enriched in Chemokine pathway. Five of the top 10 hub genes (AKT1, TNF, MAPK3, JUN, MAPK1) are enriched in T cell receptor pathway. Four of the top 10 hub genes (AKT1, MAPK3, JUN, MAPK1) are enriched in B cell receptor pathway. The three pathways are not enriched including Endocrine and other factor-regulated calcium reabsorption, Mineral absorption, and Wnt pathway.

To search for potentially relevant pathways, we obtained 50 enriched pathways as osteoporosis-related pathways and used Cytoscape software to construct the “pathway-gene” network. Studies have shown that there are many pathways associated with osteoporosis, which are mainly related to bone metabolism, inflammatory response, immune response, endocrine metabolism and apoptosis (Liu et al., 2020). Based on these analyses, Gusuibu may exert anti-osteoporosis effects through many pathways and hub genes.

As an osteoblast-specific transcription factor OSX is required for bone formation and osteoblast differentiation (Nakashima et al., 2002; Zhang et al., 2008). OSX was originally discovered as a BMP2 inducible gene in C2C12 mesenchymal stem cells, and OSX knock-out mice lack bone formation (Nakashima et al., 2002). OCN is a well-known downstream osteoblastic gene of OSX. There are three critical transcription factors required for osteoblast differentiation and bone formation: IHH, RUNX2 and OSX. OSX is the only osteoblast-specific transcription factor identified so far, and it is only expressed in cells in the bone matrix and the inner (endosteum) and outer (periosteum) bone surfaces. OSX discovery as a master regulator of bone formation and osteoblast differentiation opens a new window to bone biology. As far as we know, this is the first study to discover that one active compound kaempferol in Gusuibu activated the expression of osteoblast specific transcription factor OSX through JNK kinase pathway.

In summary, through network pharmacology, we explored the potential mechanisms by which Gusuibu exerts its anti-osteoporosis effects with a multi-component, multi-target, and multi-pathway profile. Potential therapeutic targets for Gusuibu to exert anti-osteoporosis effects include AKT1, IL6, JUN, TNF, MAPK3, VEGFA, EGFR, MAPK1, CASP3, PTGS2. Most of the key active compounds in Gusuibu are flavonoids, mainly including Luteolin, Naringenin, Kaempferol, and Beta-sitosterol. Interestingly, kaempferol was identified as an upregulator of OSX.