Abstract
Osteosarcoma is one of most malignant bone tumors in children, characterized by high recurrence and metastasis. Plumbagin refers to a bioactive compound that is isolated the herb plant of from Plumbagozeylanica zeylanica L., and it has been proven with potential anti-tumor benefits, including osteosarcoma. However, its pharmacological mechanism remains unclear comprehensively. Thus, this study aimed to reveal potential targets and molecular mechanisms in plumbagin for treating osteosarcoma through bioinformatics method and computational validation. In total, respective 379, 2727 and 2166 genes were ascertained as target genes of plumbagin, osteosarcoma and autophagy. A total of 40 shared genes were identified among plumbagin, osteosarcoma and autophagy. Further, additional 10 core genes were identified and used for enrichment analysis. The findings highlighted the regulatory actions of plumbagin on protein-binding, regulation of autophagy for playing anti-osteosarcoma role. Enriched pathway analysis findings disclosed main molecular pathways, including microRNAs in cancer signaling pathway, Notch signaling pathway. Molecular docking data found that the optimal docking affinity and binding energy between plumbagin and scored protein receptors of glycogen synthase kinase 3 beta (GSK3B), histone deacetylase 2 (HDAC2), poly (ADP-ribose) polymerase 1 (PARP1). Our preclinical study investigates the possible therapeutic mechanism of plumbagin against osteosarcoma, indicating that plumbagin exhibited anti-osteosarcoma features via regulation of core target genes associated with autophagy. Current research findings may provide the scientific ideas and evidences for screening bioactive compound against osteosarcoma.
Introduction
Osteosarcoma represents an infrequent pediatric malignancy that will develop to be lethal outcome if untreated medically. 1 In the past 30 years, the event-free survival rate of osteosarcoma is still not achieved notably and its burden is heavy to patient or society. 2 Management of osteosarcoma seems to be challenges as its histopathologic patterns are diverse and indicative biomarkers are deficient for most cases. 3 Osteosarcoma can be diagnosed globally, and it can be more prevalent in less-developed countries with increasing cases. 4 The incidence of osteosarcoma is found to be regional epidemiologically, and susceptibility of gene-related osteosarcoma is characteristic. 5 Osteosarcoma should be diagnosed accurately before being effectively treated by clinical standards. 6 Medical diagnosis of osteosarcoma is timely required owing to late disease identification, such as using plain radiography. 7 Encouragingly, adjuvant chemo-treatment combined with neoplastic resection can contribute to enhanced overall 10-year survival rate from 30% to around 50% in patients with osteosarcomas. 8 Clinical chemotherapy treating osteosarcomas is commonly prescribed, including cisplatin and ifosfamide. 9 However, drug resistance or reduced therapeutic effectiveness progressively occurs in long-range chemotherapy against osteosarcoma. 10 Hence, the oncologist should be actively encouraged to dig out innovative strategy treating osteosarcoma.
Natural compounds isolated from medicinal herbs are gained increasing traction as parts of bioactive compounds are preclinically found with pharmacological activities. Plumbagozeylanica zeylanica L is a traditional herb that found with promising pharmacological actions. 11 Experimental data report that extracts from Plumbagozeylanica zeylanica L play the antibacterial effect via in vitro assessment. 12 Plumbagin refers to one of main phytochemicals extracted from Plumbagozeylanica zeylanica L. Pharmacologically, plumbagin is preclinically proven to have potent anti-proliferative effects against human cancers. 13 Although the limited solubility and bioavailability of plumbagin are found, the optimal preparation formulation or structural modification will be developed for effective drug to its clinical application against cancer. 14 It is reported that improved intestinal absorption and solubility of plumbagin via lipid-based excipient Gelucire for increasing the extraction efficiency and bioavailability. 15 In underlying mechanisms, anti-cancer action of plumbagin may be related to firsthand induction of apoptosis, autophagy, cell cycle arrest, reactive oxygen species (ROS)-cytotoxicity and exertion of anti-angiogenic pathways. 16 Previous cell culture studies show that plumbagin elicits the anti-neoplastic action against human osteosarcoma cells through inducing cancer cell death.17,18 Currently, comprehensive assessment of action target and molecular mechanism of plumbagin treating osteosarcoma is still unclarified. Network pharmacology-based bioinformatic analysis refers to an efficient approach to investigate main target and mechanisms of candidate drugs treating of diseases, as highlighted in previous studies.19,20 Current research employed network pharmacology and molecular docking methods to identify the anti-osteosarcoma targets and mechanisms related to plumbagin action via regulating autophagy. Then, the network pharmacology findings were further assessed through molecular docking validation to exhibit candidate therapy strategy against osteosarcoma.
Materials and Methods
Potential Drug Target Genes for Plumbagin
The molecular structure of plumbagin (C11H8O3) was downloaded from PubChem database (https://pubchem.ncbi.nlm.nih.gov). Representational drug target genes of plumbagin were discerned from accessible databases of SuperPred (http://bioinformatics.charite.de/superpred), Pharm-Mapper (http://lilab.ecust.edu.cn/pharmmapper) and Swiss Target Prediction (http://www.swisstargetprediction.ch/index.php). These common targets attained from dissimilar databases were united to produce the target genes of plumbagin against osteosarcoma.
Potential Disease Target Genes for Osteosarcoma and Autophagy
Potential target genes of osteosarcoma and autophagy were gotten through searching available databases of GeneCards (https://www.genecards.org/). Venn diagram analysis was used to analyze overlapped target genes and identify the potential genes of plumbagin against osteosarcoma associated with autophagy.
Generation of Protein-Protein Interaction (PPI) Network Map
In methodology, PPI assay is a promising strategy to comprehensively determine biological functions and interaction network among overlapped proteins/genes before eventually being screened out core targets. These potential targets of plumbagin against osteosarcoma were imported to String database (https://string-db.org/) to produce a PPI picture, in which the “Organism” was applied with “Homo sapiens,” and the minimum interaction score was set to standardize medium confidence with 0.400. Then, the data using tsv format were submitted to Cytoscape software (version 3.8.2) for functional enrichment visualization and identification of gene-gene relations (core target genes).
Functional Annotation and Enrichment Analyses
All identified core target genes were remitted to the Database for Annotation, Visualization and Integrated Discovery database (http://david.abcc.ncifcrf.gov/) for performing Gene Ontology (GO)-functional annotation and Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway enrichment analyses, respectively. These KEGG entries were identified with the significance scores of P < .01 for assay, before being classified via fold algorithm in enrichments. As a result, GO annotations were grouped as following categories: biological process (BP), cellular component (CC), and molecular function (MF). And main GO terms and KEGG pathways were optimized and visualized by using P-value (P < .01), respectively. To better reveal anti-osteosarcoma target and mechanism of plumbagin, a comprehensive network was constructed by applying Cytoscape software (version 3.8.2).
Molecular Docking Assay of Core Target Proteins
Methodologically, the target proteins (core genes) and plumbagin were applied as receptors and ligand respectively. The molecular structure file of plumbagin (ligand) was downloaded from PubChem database, and the structures of GSK3B, HDAC2, PARP1 proteins (receptors) were gotten from Protein Data Bank (PDB; http://www.rcsb.org/pdb/) database. And then these data were preprocessed with PyMol software (https://www.pymol.org/), a molecular visualization system. Then, AutoDock software (https://autodock.scripps.edu/) was used to directionally add polar hydrogens, charges to target receptors before the files of ligand and receptor being converted as pdbqt format. And the active/functional pockets within target receptors were identified and conserved as gpf files. Additionally, AutoDock vina programs were used to highlight molecular docking model features, and the docking data were annotated and visualized through using PyMol software. The binding energy was determined for binding competency between ligand and receptors. It is noted that a binding energy less than −5 kcal/mol suggests the stable and effective affinity capability among ligand and receptors.
Results
Preparation of Candidate Target Genes
After screening by using online database platform, all 378 plumbagin-related target genes were gotten (Figure 1A). And other 2727 osteosarcoma-related and 2166 autophagy-related genes were acquired after being removed duplicate targets (Figure 1A). These cluster target genes identified were overlapped among plumbagin, osteosarcoma, autophagy, and then total 40 shared targets were discovered in a Venn diagram highlighting in relational network map (Figure 1B).

(A) Venn diagram of candidate therapeutic target genes of plumbagin against osteosarcoma relating to autophagy. (B) Network diagram of protein-protein interactions among shared genes. (C) Screening and identifying of core genes before construction of an interaction network.
Generation of Specialized PPI Network
All 40 shared targets were submitted to the String database (http://string-db.org) to produce a PPI network map. Then, top target genes with greatest degree scores were determined as core targets, characterized with correlated individually in network map. These core targets that potentially acted on plumbagin against osteosarcoma included TP53, glyceraldehyde 3-phosphate dehydrogenase (GAPDH), transcription 3 (STAT3), heat shock protein 90 alpha family class A member 1 (HSP90AA1), E1A binding protein P300 (EP300), HDAC1, HDAC2, PARP1, peroxisome proliferator-activated receptor gamma (PPARG), GSK3B, histone deacetylase 6 (HDAC6), B-cell lymphoma-2 like 1 (BCL2L1), HDAC4, Notch1 (NOTCH1), cyclin-A2 (CCNA2), HDAC3, Aurora kinase A (AURKA), HDAC5, retinoic acid receptor alpha (RARA), heme-oxygenase 1 (HMOX1), cyclin-dependent kinase 5 (CDK5), respectively (Figure 1C).
Revealment of the Biological Functions in Core Target Genes
All these core genes were applied for GO and KEGG enrichment assays to systematically reveal the biological functions, in which GO-enriched terms were assigned as following parts, including BP, MF, and CC. The BP part was largely linked to BP of rhythmic process, regulation of mitochondrion organization, regulation of binding, regulation of autophagy, protein localization to nucleus, protein deacylation, protein deacetylation, macromolecule deacylation, histone deacetylation, circadian rhythm (Figure 2A). These findings certified the regulatory action of plumbagin on protein-binding functions, regulation of autophagy. The MF part was mostly related to histone deacetylase binding, RNA polymerase II-specific DNA-binding transcription factor binding, histone deacetylase activity, protein lysine deacetylase activity, DNA-binding transcription factor binding, deacetylase activity, transcription coregulator binding, hydrolase activity, acting on carbon-nitrogen (but not peptide) bonds, in linear amides, NF-kappaB binding, tau protein binding (Figure 2B). The CC part was mainly related to histone deacetylase complex, transcription repressor complex, pronucleus, NuRD complex, CHD-type complex, neuronal cell body, Sin3 complex, Sin3-type complex, cyclin-dependent protein kinase holoenzyme complex, replication fork (Figure 2C). And then, KEGG-enriched analysis findings were correlation of Viral carcinogenesis, Thyroid hormone signaling pathway, MicroRNAs in cancer, Transcriptional misregulation in cancer, Homo sapiens (human), Alcoholism, Neutrophil extracellular trap formation, Human papillomavirus infection, Lipid and atherosclerosis, Notch signaling pathway (Figure 2D). The pharmacological role of plumbagin contributed to regulation of microRNA metabolism and Notch signaling pathway for potentially treating cancers. Together, these bioinformatics findings revealed the regulator mechanism of action of plumbagin against osteosarcoma.

GO-based enrichment analysis findings of core genes, as exhibited in BP (A), MF (B), CC (C). KEGG-based molecular pathways showing the core targets in plumbagin against osteosarcoma (P < .05) (D).
Verification With Molecular Docking Data
For the purpose of validating the stably binding feature between plumbagin (ligand) and GSK3B, HDAC2, PARP1 target proteins (receptors), molecular docking analysis then was performed. According to threshold values of effectively docking, binding energy less than −5.0 kcal/mol suggested that the ligand was well docked with the target receptor. The binding energy value between plumbagin and 1PYX in GSK3B was −7 kcal/mol, shaping hydrogen bonds at amino acid residues of ASP-133, VAL-135 (Figure 3A). Plumbagin docked to 3MAX in HDAC2 with a binding energy value was −8.2 kcal/mol, building hydrogen bond at amino acid residue of ALA-124 (Figure 3B). The binding energy value between plumbagin and 3GJW in PARP1 was −7.7 kcal/mol, forming hydrogen bond at amino acid residue of GLY-202 (Figure 3C).

Visualization of molecular docking feature and binding position of plumbagin with these core target proteins, including plumbagin-1PYX (GSK3B; A), plumbagin-3MAX (HDAC2; B), and plumbagin-3GJW (PARP1; C).
Discussion
Osteosarcoma refers to one of most frequent malignancies of primary bone sarcomas, characterized by a high proclivity for potent invasiveness and uncontrolled metastasis. 21 Although surgical intervention united with pharmacotherapy results in vastly meliorated the efficacy against osteosarcoma, the elimination of metastasis and recurrence of osteosarcomas remains unsatisfying. 22 Moreover, the effectiveness of existing drug therapy occurs with unavoidable medical challenge as targeted treatment and effective bioavailability are insufficient, and poor susceptibility and progressive drug resistance are also reported clinically. 23 The growing studies show that cancer cell autophagy may affect tumor growth and malignant progression. 24 Other increasing evidences suggest that autophagy may be a potential anti-osteosarcoma target for screening candidate agent, as autophagy may relate to the roles of regulating cell proliferation, metastasis, chemotherapy and radiotherapy. 25 Some of natural bioactive ingredients, such as quinacrine, are found with anti-proliferative effect against osteosarcoma through targeting of autophagy. 26 In this study, we aimed to explore the anti-osteosarcomas target and mechanism relating to autophagy of plumbagin through using the network pharmacology and molecular docking analyses. Following with network pharmacology analysis, total 40 shared target genes were gotten via merging and optimizing the genes of plumbagin and osteosarcoma relating to autophagy. Furthermore, TP53, GAPDH, STAT3, HSP90AA1, EP300, HDAC1, HDAC2, PARP1, PPARG, GSK3B, HDAC6, BCL2L1, HDAC4, NOTCH1, CCNA2, HDAC3, AURKA, HDAC5, RARA, HMOX1, CDK5 were identified as core target genes before being performed with functional enrichment analysis and molecular docking validation. The KEGG analysis findings indicated that top molecular pathways may include microRNAs in cancer signaling pathway, Notch signaling pathway. It is found that certain endogenous MicroRNAs (miRNAs) are differentially expressed in human osteosarcoma or cancer cell line samples, characterized with malignant phenotype and adverse prognosis. Targeting of these small-molecule miRNAs and related signaling pathway can be used for exploring and developing anti-osteosarcoma agent exploring. 27 It is reported that dysregulated miRNAs can be involved in induction of cancerous development from normal cells through canonical autophagy pathway, and subsequently cause chemoresitance. 28 The Notch signaling pathway exerts a crucial role in the occurrence and progression of human cancer. Clinical observations report that activated Notch pathway is detected in human osteosarcoma samples and correlated to high recurrence and poor prognosis. Thus, targeting of Notch signaling-based treatment may be one of anti-osteosarcoma strategies. 29 Autophagy-related gene is elevated expression in osteosarcoma samples, and induces cell epithelial-to-mesenchymal transition via activating the Notch signaling pathway. 30 The molecular docking validation findings revealed that the binding affinity and docking energy between plumbagin (ligand) and core target receptors were effective and stabilized, including GSK3B, HDAC2, PARP1 target proteins. Increasing research data indicate that GSK3B (glycogen synthase kinase 3β) mediates a key role in oncogenesis. Thus, it is concluded that GSK3B inhibitor can be used for target treatment against malignant cancers. 31 GSK3B-medicated tumorigenesis through phosphorylation of Unc-51-like autophagy activating kinase 1 (ULK1) is involved in induction of autophagy. 32 GSK3B is highly expressed in human osteosarcoma samples, and then targeting GSK3B therapy contributes to playing anti-osteosarcoma effectiveness followed with minimal adverse effects. 33 The growing data indicate that HDAC2 may serve as a promising prognostic biomarker for screening malignant tumors, and then HDAC2 may be as a potential target for developing anti-cancer agent. 34 Abnormal expression of HDAC2 plays an important role in neoplastic development via driving autophagy-related genes. Therefore, targeting of HDAC2 may be used for prophylaxis and treatment of malignant cancers. 35 Highly expressed HDAC2 is found in human osteosarcoma samples, and resulted in reduced overall survival. And targeting of HDAC2 inhibition is reported with reduction of osteosarcoma cell invasion and migration. 36 PARP1 (Poly(ADP-ribose) polymerase 1) is pivotal for the occurrence and development of human malignant cancers. 37 PARP1-mediated autophagy is a key gateway for drug resistance when using cancer pharmacotherapy. 38 It is reported therapeutic vulnerability of PARP1 inhibitor in RB1-mutant sporadic osteosarcoma. 39 To sum up, integrated bioinformatics- and reference-based data speculated that GSK3B, HDAC2, PARP1 may be therapeutic targets of plumbagin against osteosarcoma relating to autophagy.
In this study, we attempted to explore the anti-osteosarcomas target and mechanism of plumbagin through preclinical data from network pharmacology and molecular docking analyses. Meanwhile, the limitations and proposed the prospect were also showed in current study. This research was designed to use bioinformatics analysis, and these computational findings should be further verified by in vitro and in vivo experiments. Moreover, there was marked diversity between clinical osteosarcoma subtypes. Notably, this study was just based on online-database analysis, and there was unable to conduct detailed analysis of different osteosarcomas subtypes. Thus, updated online-databases will be further determined in future before experimental validation.
Conclusion
Taken together, this research aims to reveal the multi-targets and multi-molecular pathways of plumbagin against osteosarcoma through using network pharmacology and molecular docking approaches. And then, the computational validation findings suggest that plumbagin may exert anti-osteosarcoma benefits through targeting regulation of GSK3B, HDAC2, PARP1 core proteins relating to autophagy. These bioinformatics results highlight that plumbagin may be a candidate anti-osteosarcoma compound before further experimental validation.
Footnotes
Author Contributions
All authors agreed to be accountable for the aspects of this research work. Muhua Liang, Rubiao Qiu proposed this study, prepared and wrote the initial draft of the manuscript, and given suggestions for revision. Rubiao Qiu, Xueyu Li, Yanjuan Li, Fufeng Yuan, Zhongxi Cen, Guoshu Huang, Xiong Chen, Chaohui Fan collected data and analyzed data, summarized the data and plotted diagrams.
Funding
The authors disclosed receipt of the following financial support for the research, authorship, and/or publication of this article: This research was supported by the Natural Science Foundation of Guangxi Province (No. 2023GXNSFAA026111), Guangxi Key Laboratory of Basic Research on Birth Defects Prevention and Control Open Subject (No. GXWCH-ZDKF-2022-24).
Declaration of conflicting interests
The authors declared no potential conflicts of interest with respect to the research, authorship, and/or publication of this article.
