Next Article in Journal
EAE of Mice: Enzymatic Cross Site-Specific Hydrolysis of H2A Histone by IgGs against H2A, H1, H2B, H3, and H4 Histones and Myelin Basic Protein
Previous Article in Journal
Structure-Based Design of Potent Peptidomimetic Inhibitors Covalently Targeting SARS-CoV-2 Papain-like Protease
Previous Article in Special Issue
Assessing the Differential Methylation Analysis Quality for Microarray and NGS Platforms
 
 
Search for Articles:
Title / Keyword
Author / Affiliation / Email
Journal
Article Type
 
 
Advanced Search
 
Section
Special Issue
Volume
Issue
Number
Page
 
Journals
IJMS
Volume 24
Issue 10
10.3390/ijms24108635
ijms-logo
Submit to this Journal Review for this Journal Propose a Special Issue
► Article Menu

Article Menu

share Share announcement Help format_quote Cite question_answer Discuss in SciProfiles thumb_up
...
Endorse textsms
...
Comment
first_page
settings
Order Article Reprints
Font Type:
Arial Georgia Verdana
Font Size:
Aa Aa Aa
Line Spacing:
Column Width:
Background:
Communication

Identification of Potential Biomarkers and Small Molecule Drugs for Bisphosphonate-Related Osteonecrosis of the Jaw (BRONJ): An Integrated Bioinformatics Study Using Big Data

by
Kumarendran Balachandran
1,
Roszalina Ramli
2,
Saiful Anuar Karsani
3 and
Mariati Abdul Rahman
1,*
1
Department of Craniofacial Diagnostics and Biosciences, Faculty of Dentistry, University Kebangsaan Malaysia, Jalan Raja Muda Abdul Aziz, Kuala Lumpur 50300, Malaysia
2
Department of Oral and Maxillofacial Surgery, Faculty of Dentistry, Universiti Kebangsaan Malaysia, Kuala Lumpur 50300, Malaysia
3
Institute of Biological Sciences, Faculty of Science, Universiti Malaya, Kuala Lumpur 50603, Malaysia
*
Author to whom correspondence should be addressed.
Int. J. Mol. Sci. 2023, 24(10), 8635; https://doi.org/10.3390/ijms24108635
Submission received: 28 February 2023 / Revised: 18 April 2023 / Accepted: 25 April 2023 / Published: 11 May 2023
(This article belongs to the Special Issue The Development and Application of Computational Models for Identifying Disease Markers)
Browse Figures
Versions Notes

Abstract

:
This study aimed to identify potential molecular mechanisms and therapeutic targets for bisphosphonate-related osteonecrosis of the jaw (BRONJ), a rare but serious side effect of bisphosphonate therapy. This study analyzed a microarray dataset (GSE7116) of multiple myeloma patients with BRONJ (n = 11) and controls (n = 10), and performed gene ontology, a pathway enrichment analysis, and a protein–protein interaction network analysis. A total of 1481 differentially expressed genes were identified, including 381 upregulated and 1100 downregulated genes, with enriched functions and pathways related to apoptosis, RNA splicing, signaling pathways, and lipid metabolism. Seven hub genes (FN1, TNF, JUN, STAT3, ACTB, GAPDH, and PTPRC) were also identified using the cytoHubba plugin in Cytoscape. This study further screened small-molecule drugs using CMap and verified the results using molecular docking methods. This study identified 3-(5-(4-(Cyclopentyloxy)-2-hydroxybenzoyl)-2-((3-hydroxybenzo[d]isoxazol-6-yl) methoxy) phenyl) propanoic acid as a potential drug treatment and prognostic marker for BRONJ. The findings of this study provide reliable molecular insight for biomarker validation and potential drug development for the screening, diagnosis, and treatment of BRONJ. Further research is needed to validate these findings and develop an effective biomarker for BRONJ.

1. Introduction

Bisphosphonate-related osteonecrosis of the jaw (BRONJ) is characterized by an exposed refractory bone in the area of the jaw, caused by an adverse effect of bisphosphonate (BP) treatment [1]. Currently, there are no effective biomarkers for BRONJ. The understanding of BRONJ pathophysiology remains theoretical. It is believed that aberrations in bone metabolism are caused by the deposition of bisphosphonates (BP) in the bone cause the homeostatic imbalance of osteoblasts and osteoclasts [2]. However, many studies have suggested a genetic link of the condition as not everyone treated with BP develops BRONJ [3,4].
Diagnostic bioinformatics is the application of bioinformatic techniques towards the development and improvement of diagnostic tests [5]. The impact of bioinformatics has been significant, as it allows for the efficient and accurate identification of diseases and the development of personalized treatment plans [6,7]. For example, cancer diagnostics using a large amount of genetic and molecular data has been used to identify genetic markers for different types of cancers [8,9]. Additionally, non-invasive prenatal testing (NIPT) has allowed for detection of chromosome abnormalities in fetal DNA [10]. In analyzing large data, bioinformatics has aided in the identification of new drug targets and the optimization of drug designs for the development of personalized medicine [11,12].
A broad spectrum of gene expression data is available online. If thoroughly investigated, these data may provide predictive insights into the molecular mechanisms underlying diseases. As such, this study aimed to utilize these data, together with powerful informatic tools, to unravel and identify potential biomarkers for BRONJ diagnosis, and to identify potential drugs and drug targets for the treatment of BRONJ.
A robust analytical method was used in this study. While an analysis of a single gene chip, as performed in this study, may yield false-positive results, we used a robust multichip array average (RMA) method, which integrated four important features: strong robustness to noise, incomplete ranking, significant scores in the result ranking, and high computational ranking. To the best of our knowledge, no previous study has used the RMA method that facilitated this study to identify differentially expressed genes (DEG) in BRONJ. A group from Korea used similar datasets for a combined predictor analysis and did not conduct any pre-processing of the data [13]. Using a non-parametric statistical analysis, they found 200 significant genes [13]. A more recent study used GSE7116 for data mining for anti-cancer therapeutics using the DGIdb database [14]. The GEO2R tool used in the previous anti-cancer therapeutic study poses several limitation [14,15]. GEO2R uses statistical tests to identify differentially expressed genes; however, these tests make certain assumptions about the data, such as the normality and homogeneity of variance, without the ability to pre-process the data. We have provided a summary of studies that have previously used GSE7116 in Table 1.
In this study, we performed gene ontology (GO), Kyoto Encyclopedia of Genes and Genomes (KEGG) enrichment, and protein–protein network analyses. We also provided the relationship between different genes and their regulatory networks, using the cytoHubba plugin of Cytoscape. Previous studies have successfully utilized this approach for biomarker discovery [20,21,22]. Furthermore, from the results of the above analysis, we sourced small drug molecules that interacted with key hub genes. To validate this, we used a molecular docking method to confirm the possible mechanism of action. The steps of this procedure are explained in our Section 4.1.

2. Results

2.1. Identification of DEGs in BRONJ

Our search strategy yielded several relevant studies, and we downloaded and analyzed the microarray dataset GSE7116. This study was a platform based on GPL570, containing total samples from 11 MM patients with BRONJ and 10 MM patients without BRONJ. A volcano plot of the microarray results after RMA normalization (p < 0.05) is shown in Figure 1. A total of 1481 DEGs (381 upregulated and 1100 downregulated) were identified.

2.2. Functional and Network Analysis of DEGs

We used the Database for Annotation, Visualization, and Integrated Discovery (DAVID) to perform a functional analysis of these DEGs. As shown in Figure 2, we found that DEGs between the control and BRONJ groups were significantly enriched in the apoptotic process, T cell receptor signaling pathway, RNA splicing, RNA binding, PI3k-Akt signaling pathway, regulation of actin cytoskeleton, MAPK signaling, and lipid and atherosclerosis pathways.

2.3. PPI Network and Hub Genes

We used STRING to examine and identify possible PPI networks. A PPI network was constructed using the identified DEGs, which is presented in Figure 3. The cytoHubba plug-in in Cytoscape was used to determine the most significant hub gene(s) from the PPI network. We expanded the Venn diagram (Figure 3F) to intersect the most significant genes according to the density of the maximum neighborhood component (MNC), maximal clique centrality (MCC), density of maximum neighborhood component (DMNC), closeness, betweenness, and degree. This led to the identification of seven hub genes (Figure 4). Functional analyses of the seven hub genes were performed using DAVID. The expression patterns of these genes are provided in the Supplementary Materials.

2.4. Small Molecule Drug Screening

We utilized the CMap network to analyze the seven identified DEGs. The top six compounds with the highest negative enrichment scores were identified as potential therapeutic targets for BRONJ, Alantolactone, 1,4-dihydronicotinamide adenine dinucleotide, Glycitein, Napabucasin, 3-(5-(4-(cyclopentyloxy)-2-hydroxybenzoyl)-2-((3-hydroxybenzo[d]isoxazol-6-yl) meth-oxy) phenyl) propanoic acid, and Tanshinone IIA. The chemical structures of these compounds, and the hydrogen bonds surrounding the carbon compound, are illustrated in Figure 5.

2.5. Verification by Molecular Docking

Using PyRx 0.8, the screened small-molecule drugs were docked with seven core targets (ACTB, STAT3, FN1, JUN, PTPRC, GAPDH, and TNF). A binding energy lower than 0 suggests that the ligand and receptor bind spontaneously. The possibility of a higher activity was determined by the stability and binding energy of the conformation. The lowest binding energy was usually less than −11 kcal mole−1, showing that the target protein had a high affinity for the active component. The summary of binding energies is shown in Figure 6. A small molecule therapeutic docking target was developed. The bonding of a drug molecule and its interaction with identified genes is shown in Figure 7. For an example, Tanshionone IIA exerts its biological activity by attaching to, and creating, hydrogen bonds with four amino acids near the active site, as represented by the dotted line.

3. Discussion

Traditional gene expression analysis poses a challenge considering the biological heterogeneity and technical biases of sequencing platforms [23]. In this study, we utilized methods to eliminate these limitations through a robust normalization and scaling based on the relative ranks of gene expression levels, as proposed by several previous studies [24,25]. In this study, we identified 1481 DEGs, including 381 upregulated and 1100 downregulated genes. Further analysis narrowed down our search to seven hub genes: FN1, TNF, JUN, STAT3, ACTB, GAPDH, and PTPRC. Through CMap and molecular docking, we found that 3-(5-(4-(cyclopentyloxy)-2-hydroxybenzoyl)-2-((3-hydroxybenzo[d]isoxazol-6-yl) methoxy) phenyl) propanoic acid (PubChem ID: 23626877) could be used to target all of the postulated genes to reverse BRONJ.
All hub genes were directly or indirectly involved in bone metabolism. The FN1 gene provides instructions for producing fibronectin-1 protein [26]. FN1 has been shown to mediate chondrocyte adhesion [27]. Upregulation of FN1 mediates fracture healing by activating the TGF-B/PI3K/Akt signaling pathway [28]. TNF is a protein-coding gene that is associated with malaria and asthma [29]. Additionally, TNF plays an important role in skeletal system-induced inflammatory processes [30]. Jun encodes proteins via GTPase activator activity. In several studies, Jun has been postulated to cause bone development [31,32,33]. ACTB provides instructions for the expression of B-actin. In a study on postmenopausal osteoporosis subjects, it was shown that ACTB was aberrated by being an unsuitable house-keeping gene for expression assays [34]. GAPDH is a moonlighting protein that serves as a glycolytic enzyme and an uracil DNA glycosylase [35]. In a study conducted to investigate the possible effects of BP on GAPDH mRNA expression, GAPDH expression decreased, similar to our findings [36]. PTPRC is a tyrosine phosphatase (PTP) protein. PTPRC is also known as the CD45 antigen or leukocyte common antigen (LCA) [37].
Based on our scoping review of the genes studied in BRONJ, STAT3 is postulated to be activated subsequent to the immune response. Treatment with BP elevates IL-6 and IL-36a levels, which then activates the STAT3 pathway [38]. The molecular function of this pathway is illustrated in Figure 8. Several studies have shown the differential pathways of nitrogenous versus non-nitrogenous BP on STAT3. Disruption of this pathway may effectively control RANKL expression.
3-{5-[4-(cyclopentyloxy)-2-hydroxybenzoyl]-2-[(3-hydroxy-1,2-benzoxazol-6-yl) methoxy] phenyl} propanoic acid is an organic compound that is classified as a benzophenone, which has a ketone attached to two phenyl groups. Despite being a relatively unknown compound, it has been detected in human blood, indicating that it is not a naturally occurring metabolite, and is only present in people who have been exposed to it or its derivatives [39]. This compound is considered part of an individual’s exposome, which encompasses all the environmental and occupational exposures that affect their health throughout their lifetime, starting before birth [39]. This study provides a basis for the further exploration and validation of 3-{5-[4-(cyclopentyloxy)-2-hydroxybenzoyl]-2-[(3-hydroxy-1,2-benzoxazol-6-yl) methoxy] phenyl} propanoic acid for therapeutic purposes.
This study lends itself to several limitations. Although bioinformatics data mining can provide valuable information, available datasets, particularly those related to BRONJ, are limited. As such, we utilized a single center study. Additionally, a single chip analysis should be interpreted with caution as it may account for limited coverage, technical variability, and cross hybridization, making it less feasible to arrive at a plausible conclusion.

4. Materials and Methods

4.1. Flowchart of Study

Figure 9 illustrates the pipeline of this study. Firstly, the raw data from a gene expression experiment or a protein study was collected and compiled in a geo dataset. Then, the data was preprocessed and analyzed using statistical methods, such as RMA analysis, to normalize and filter the data. Furthermore, gene ontology (GO) and a Kyoto Encyclopedia of Genes and Genomes (KEGG) analysis were performed to identify the biological processes, molecular functions, and pathways that were significantly enriched in the dataset. Subsequently, a protein–protein interaction (PPI) network was constructed to visualize the interactions between different proteins and to identify the highly connected hub genes that play a crucial role in the biological system being studied. Once the hub genes were identified, small drug molecules were screened and identified using virtual screening methods to target the protein products of these hub genes. A molecular docking verification was then performed to predict the binding affinity and the stability of the drug–protein complex.

4.2. Differential Expression Analysis

The expression profile for BRONJ was obtained from the Gene Expression Omnibus (GEO) database [(https://www.ncbi.nlm.nih.gov/, accessed on 28 November 2022)]. The search strategy [“bisphosphonate related osteonecrosis of the jaw” (MeSH Terms) OR “bisphosphonate related osteonecrosis of the jaw” (All fields)] and [“Home sapiens” (Organism) and “Expression profiling” (Filter)] were used. The data contained 11 multiple myeloma (MM) patients with BRONJ, against 10 age-matched MM patients without BRONJ. Further details are provided in Table 2. The ‘limma’ package was used to normalize the data in R using the RMA method (version 4.2.2, http://www.R-project.org/, accessed on 28 November 2022). Using the R package ‘limma,’ a linear model was used to assess differential expression. The cut-off parameters for determining DEGs were [log2 fold change (FC)] > 1 and p-value p < 0.05.

4.3. Gene Ontology and Pathway Analysis

GO enrichment and KEGG pathway analyses were performed using a freely available Database for Annotation, Visualization, and Integrated Discovery (DAVID) (DAVID 6.8, http://david.ncifcrf.gov/, accessed on 28 November 2022). GO annotates gene products and functions into three categories: biological processes (BP), molecular functions (MF), and cellular components (CC), allowing us to compare and analyze genes across control versus BRONJ subjects. A p-value of 0.05 was used as the cutoff criterion.

4.4. Protein–Protein Interaction (PPI) Network Analysis

The STRING online database was used to construct a functional PPI network for the identified DEGs (https://string-db.org/, accessed on 28 November 2022). The cut-off point was set to a credibility score higher than 0.4. PPI visualization was performed using the Cytoscape software.

4.5. Hub Gene Mapping

The Cytoscape cytoHubba plugin was used to select hubs from the DEGs. (https://apps.cytoscape.org/apps/cytohubba, accessed on 28 November 2022). Cytoscape is an open-source software platform used for visualizing and analyzing molecular interaction networks. CytoHubba is a plugin for Cytoscape that provides various network analysis algorithms to identify important nodes or subnetworks within a given network. Several topological algorithms were used to identify interactome networks in this plug-in. The following metrics were used: density of the maximum neighborhood component (MNC), maximal clique centrality (MCC), density of the maximum neighborhood component (DMNC), closeness, betweenness, and degree. Venn diagrams were used to locate the shared genes.

4.6. Identification of Small Molecule drug Candidates

To assess the potential effectiveness of the medication treatment for BRONJ, we measured the binding energy between the small-molecule medications predicted in CMap and the possible target proteins of BRONJ. The crystal structures of the principle targets were obtained from the RCSB Protein Data Bank (PDB, http://www.rcsb.org, accessed on 28 November 2022), whereas the mol2 file formats of the compounds were retrieved from the PubChem database. Target proteins were extracted from PyMOL 2.3.2, and their ligands were stored in PDB format. The target protein was then hydrogenated, and its charge was calculated and stored in PDBQT format using AutoDock Tools 1.5.6 software. Grid box data for the protein of interest were obtained, and molecular docking was performed using AutoDock Vina 1.1.2. The results of molecular docking were visualized using PyRx software. PyRx 0.8 is a useful tool for computational drug discovery and design, used to perform virtual screening experiments by docking small molecule ligands to protein targets.

5. Conclusions

There is currently no known biomarker for BRONJ. The present study used bioinformatics to analyze the gene chip data of MM patients with BRONJ. The resulting findings allowed us to identify potential molecular biomarkers for BRONJ. Apart from that, we were also able to identify potential therapeutic targets for BRONJ. All identified hub genes were either directly or indirectly related to bone metabolism, which could potentially be aberrated in disease, and hence a factor in the development of BRONJ. The present results are based only on bioinformatics analysis, and further validation is required to verify this finding. Further efforts in utilizing microarray studies may benefit efficient biomarker development.

Supplementary Materials

The supporting information can be downloaded at: https://www.mdpi.com/article/10.3390/ijms24108635/s1.

Author Contributions

Conceptualization, K.B. and M.A.R.; methodology, K.B.; software, K.B.; validation, K.B.; formal analysis, K.B.; writing—original draft preparation, K.B.; writing—review and editing, K.B., M.A.R., S.A.K. and R.R.; supervision, M.A.R., S.A.K. and R.R. All authors have read and agreed to the published version of the manuscript.

Funding

This research was funded by the Malaysia Ministry of Higher Education through the Fundamental Research Grants Scheme (FRGS), grant number FRGS/1/2020/SKK0/UKM/02/27.

Institutional Review Board Statement

Not applicable.

Informed Consent Statement

Not applicable.

Data Availability Statement

All data can be easily obtained and linked in the respective sections.

Conflicts of Interest

The authors declare no conflict of interest. The funders had no role in the design of this study; in the collection, analyses, or interpretation of data; in the writing of the manuscript; or in the decision to publish the results.

References

  1. Ruggiero, S.L.; Dodson, T.B.; Aghaloo, T.; Carlson, E.R.; Ward, B.B.; Kademani, D. American Association of Oral and Maxillofacial Surgeons’ Position Paper on Medication-Related Osteonecrosis of the Jaws—2022 Update. J. Oral Maxillofac. Surg. 2022, 80, 920–943. [Google Scholar] [CrossRef] [PubMed]
  2. Allen, M.R.; Burr, D.B. The Pathogenesis of Bisphosphonate-Related Osteonecrosis of the Jaw: So Many Hypotheses, So Few Data. J. Oral Maxillofac. Surg. 2009, 67, 61–70. [Google Scholar] [CrossRef] [PubMed]
  3. Stains, J.P.; Watkins, M.P.; Grimston, S.K.; Hebert, C.; Civitelli, R. Molecular Mechanisms of Osteoblast/Osteocyte Regulation by Connexin43. Calcif. Tissue Int. 2013, 94, 55–67. [Google Scholar] [CrossRef]
  4. Katz, J.; Gong, Y.; Salmasinia, D.; Hou, W.; Burkley, B.; Ferreira, P.; Casanova, O.; Langaee, T.Y.; Moreb, J.S. Genetic Polymorphisms and Other Risk Factors Associated with Bisphosphonate Induced Osteonecrosis of the Jaw. Int. J. Oral Maxillofac. Surg. 2011, 40, 605–611. [Google Scholar] [CrossRef] [PubMed]
  5. Wani, Y.; Rani, S.; Mir, M.R.; Sahaf, K.A.; Dar, K.A.; Ganie, N.A.; Mehraj, S.; Baqual, M.F. Advances and Applications of Bioinformatics in Various Fields of Life. Int. J. Fauna Biol. Stud. 2018, 5, 3–10. [Google Scholar]
  6. Fernald, G.H.; Capriotti, E.; Daneshjou, R.; Karczewski, K.J.; Altman, R.B. Bioinformatics Challenges for Personalized Medicine. Bioinformatics 2011, 27, 1741. [Google Scholar] [CrossRef]
  7. Kuznetsov, V.; Lee, H.K.; Maurer-Stroh, S.; Molnár, M.J.; Pongor, S.; Eisenhaber, B.; Eisenhaber, F. How Bioinformatics Influences Health Informatics: Usage of Biomolecular Sequences, Expression Profiles and Automated Microscopic Image Analyses for Clinical Needs and Public Health. Health Inf. Sci. Syst. 2013, 1, 2. [Google Scholar] [CrossRef]
  8. Sarhadi, V.K.; Armengol, G. Molecular Biomarkers in Cancer. Biomolecules 2022, 12, 1021. [Google Scholar] [CrossRef]
  9. Sokolenko, A.P.; Imyanitov, E.N. Molecular Diagnostics in Clinical Oncology. Front. Mol. Biosci. 2018, 5, 76. [Google Scholar] [CrossRef]
  10. Duboc, V.D.S.; Pratella, D.; Milanesio, M.; Boudjarane, J.; Descombes, S.D.S.; Paquis-Flucklinger, V.D.S.; Bottini, S. NiPTUNE: An Automated Pipeline for Noninvasive Prenatal Testing in an Accurate, Integrative and Flexible Framework. Brief. Bioinform. 2022, 23, bbab380. [Google Scholar] [CrossRef]
  11. Hassan, M.; Awan, F.M.; Naz, A.; Deandrés-Galiana, E.J.; Alvarez, O.; Cernea, A.; Fernández-Brillet, L.; Fernández-Martínez, J.L.; Kloczkowski, A. Innovations in Genomics and Big Data Analytics for Personalized Medicine and Health Care: A Review. Int. J. Mol. Sci. 2022, 23, 4645. [Google Scholar] [CrossRef]
  12. Wooller, S.K.; Benstead-Hume, G.; Chen, X.; Ali, Y.; Pearl, F.M.G. Bioinformatics in Translational Drug Discovery. Biosci. Rep. 2017, 37, 20160180. [Google Scholar] [CrossRef] [PubMed]
  13. Kim, K.Y.; Zhang, X.; Cha, I.H. Identifying a Combined Biomarker for Bisphosphonate-Related Osteonecrosis of the Jaw. Clin. Implant Dent. Relat. Res. 2018, 20, 191–198. [Google Scholar] [CrossRef]
  14. Zhuang, J.; Zu, J.; Zhou, C.; Sun, Y.; Kong, P.; Jing, Y. Bioinformatic Data Mining for Candidate Drugs Affecting Risk of Bisphosphonate-Related Osteonecrosis of the Jaw (BRONJ) in Cancer Patients. Dis. Markers 2022, 2022, 3348480. [Google Scholar] [CrossRef] [PubMed]
  15. Ritchie, M.E.; Phipson, B.; Wu, D.; Hu, Y.; Law, C.W.; Shi, W.; Smyth, G.K. Limma Powers Differential Expression Analyses for RNA-Sequencing and Microarray Studies. Nucleic Acids Res. 2015, 43, e47. [Google Scholar] [CrossRef] [PubMed]
  16. Leng, D.; Miao, R.; Huang, X.; Wang, Y. In Silico Analysis Identifies CRISP3 as a Potential Peripheral Blood Biomarker for Multiple Myeloma: From Data Modeling to Validation with RT-PCR. Oncol. Lett. 2018, 15, 5167–5174. [Google Scholar] [CrossRef]
  17. Sun, J.; Wen, X.; Jin, F.; Li, Y.; Hu, J.; Sun, Y. Bioinformatics Analyses of Differentially Expressed Genes Associated with Bisphosphonate-Related Osteonecrosis of the Jaw in Patients with Multiple Myeloma. Onco Targets Ther. 2015, 8, 2681–2688. [Google Scholar] [CrossRef]
  18. He, J.; Zhou, Q.; Jia, X.; Zhou, P.; Chen, L. Immune-Related Expression Profiles of Bisphosphonates-Related Osteonecrosis of the Jaw in Multiple Myeloma. Pharmazie 2021, 76, 159–164. [Google Scholar] [CrossRef]
  19. Ma, H.; Zhang, W.; Shi, J. Differentially Expressed Genes Reveal the Biomarkers and Molecular Mechanism of Osteonecrosis. J. Healthc. Eng. 2022, 2022, 8684137. [Google Scholar] [CrossRef]
  20. Zhang, X.; Wang, Z.; Hu, L.; Shen, X.; Liu, C. Identification of Potential Genetic Biomarkers and Target Genes of Peri-Implantitis Using Bioinformatics Tools. Biomed Res. Int. 2021, 2021, 1759214. [Google Scholar] [CrossRef]
  21. Zhao, C.; Quan, X.; He, J.; Zhao, R.; Zhang, Y.; Li, X.; Sun, S.; Ma, R.; Zhang, Q. Identification of Significant Gene Biomarkers of Low Back Pain Caused by Changes in the Osmotic Pressure of Nucleus Pulposus Cells. Sci. Rep. 2020, 10, 3708. [Google Scholar] [CrossRef] [PubMed]
  22. Qian, L.; Xia, Z.; Zhang, M.; Han, Q.; Hu, D.; Qi, S.; Xing, D.; Chen, Y.; Zhao, X. Integrated Bioinformatics-Based Identification of Potential Diagnostic Biomarkers Associated with Diabetic Foot Ulcer Development. J. Diabetes Res. 2021, 2021, 5445349. [Google Scholar] [CrossRef]
  23. Leek, J.T.; Scharpf, R.B.; Bravo, H.C.; Simcha, D.; Langmead, B.; Johnson, W.E.; Geman, D.; Baggerly, K.; Irizarry, R.A. Tackling the Widespread and Critical Impact of Batch Effects in High-Throughput Data. Nat. Rev. Genet. 2010, 11, 733–739. [Google Scholar] [CrossRef]
  24. Popovici, V.; Budinska, E.; Tejpar, S.; Weinrich, S.; Estrella, H.; Hodgson, G.; Van Cutsem, E.; Xie, T.; Bosman, F.T.; Roth, A.D.; et al. Identification of a Poor-Prognosis BRAF-Mutant-like Population of Patients with Colon Cancer. J. Clin. Oncol. 2012, 30, 1288–1295. [Google Scholar] [CrossRef] [PubMed]
  25. Li, B.; Cui, Y.; Diehn, M.; Li, R. Development and Validation of an Individualized Immune Prognostic Signature in Early-Stage Nonsquamous Non–Small Cell Lung Cancer. JAMA Oncol. 2017, 3, 1529. [Google Scholar] [CrossRef] [PubMed]
  26. Castelletti, F.; Donadelli, R.; Banterla, F.; Hildebrandt, F.; Zipfel, P.F.; Bresin, E.; Otto, E.; Skerka, C.; Renieri, A.; Todeschini, M.; et al. Mutations in FN1 Cause Glomerulopathy with Fibronectin Deposits. Proc. Natl. Acad. Sci. USA 2008, 105, 2538–2543. [Google Scholar] [CrossRef] [PubMed]
  27. Cadoff, E.B.; Sheffer, R.; Wientroub, S.; Ovadia, D.; Meiner, V.; Schwarzbauer, J.E. Mechanistic Insights into the Cellular Effects of a Novel FN1 Variant Associated with a Spondylometaphyseal Dysplasia. Clin. Genet. 2018, 94, 429. [Google Scholar] [CrossRef] [PubMed]
  28. Zhang, H.; Chen, X.; Xue, P.; Ma, X.; Li, J.; Zhang, J. FN1 Promotes Chondrocyte Differentiation and Collagen Production via TGF-β/PI3K/Akt Pathway in Mice with Femoral Fracture. Gene 2021, 769, 145253. [Google Scholar] [CrossRef]
  29. Osta, B.; Benedetti, G.; Miossec, P. Classical and Paradoxical Effects of TNF-α on Bone Homeostasis. Front. Immunol. 2014, 5, 48. [Google Scholar] [CrossRef]
  30. Zhao, B. TNF and Bone Remodeling. Curr. Osteoporos. Rep. 2017, 15, 126. [Google Scholar] [CrossRef]
  31. Lerbs, T.; Cui, L.; Muscat, C.; Saleem, A.; van Neste, C.; Domizi, P.; Chan, C.; Wernig, G. Expansion of Bone Precursors through Jun as a Novel Treatment for Osteoporosis-Associated Fractures. Stem Cell Rep. 2020, 14, 603–613. [Google Scholar] [CrossRef] [PubMed]
  32. Sanghani-Kerai, A.; McCreary, D.; Lancashire, H.; Osagie, L.; Coathup, M.; Blunn, G. Stem Cell Interventions for Bone Healing: Fractures and Osteoporosis. Curr. Stem Cell Res. Ther. 2018, 13, 369–377. [Google Scholar] [CrossRef] [PubMed]
  33. Haraguchi, R.; Kitazawa, R.; Imai, Y.; Kitazawa, S. Growth Plate-Derived Hedgehog-Signal-Responsive Cells Provide Skeletal Tissue Components in Growing Bone. Histochem. Cell Biol. 2018, 149, 365–373. [Google Scholar] [CrossRef] [PubMed]
  34. de Lima, C.A.D.; de Lima, S.C.; Barbosa, A.D.; Sandrin-Garcia, P.; de Barros Pita, W.; de Azevêdo Silva, J.; Crovella, S. Postmenopausal Osteoporosis Reference Genes for QPCR Expression Assays. Sci. Rep. 2019, 9, 16533. [Google Scholar] [CrossRef] [PubMed]
  35. Colell, A.; Green, D.R.; Ricci, J.E. Novel Roles for GAPDH in Cell Death and Carcinogenesis. Cell Death Differ. 2009, 16, 1573–1581. [Google Scholar] [CrossRef]
  36. Valenti, M.T.; Bertoldo, F.; Dalle Carbonare, L.; Azzarello, G.; Zenari, S.; Zanatta, M.; Balducci, E.; Vinante, O.; Lo Cascio, V. The Effect of Bisphosphonates on Gene Expression: GAPDH as a Housekeeping or a New Target Gene? BMC Cancer 2006, 6, 49. [Google Scholar] [CrossRef]
  37. Al Barashdi, M.A.; Ali, A.; McMullin, M.F.; Mills, K. Protein Tyrosine Phosphatase Receptor Type C (PTPRC or CD45). J. Clin. Pathol. 2021, 74, 548. [Google Scholar] [CrossRef]
  38. He, L.; Sun, X.; Liu, Z.; Qiu, Y.; Niu, Y. Pathogenesis and Multidisciplinary Management of Medication-Related Osteonecrosis of the Jaw. Int. J. Oral Sci. 2020, 12, 30. [Google Scholar] [CrossRef]
  39. Barupal, D.K.; Fiehn, O. Generating the Blood Exposome Database Using a Comprehensive Text Mining and Database Fusion Approach. Environ. Health Perspect. 2019, 127, 097008. [Google Scholar] [CrossRef]
  40. Raje, N.; Woo, S.B.; Hande, K.; Yap, J.T.; Richardson, P.G.; Vallet, S.; Treister, N.; Hideshima, T.; Sheehy, N.; Chhetri, S.; et al. Clinical, Radiographic, and Biochemical Characterization of Multiple Myeloma Patients with Osteonecrosis of the Jaw. Clin. Cancer Res. 2008, 14, 2387–2395. [Google Scholar] [CrossRef]
Ijms 24 08635 g001 550
Figure 1. Identification of DEGs data from GEO dataset GSE7116. Volcano plot representing differential gene expression analysis between condition A and condition B. Genes that are significantly up-regulated (log2 fold change >1, adjusted p-value <0.05) are shown in red, while genes that are significantly down-regulated (log2 fold change < −1, adjusted p-value < 0.05) are shown in blue. The x-axis represents the log2 fold change in gene expression between BRONJ vs. control, while the y-axis represents the negative logarithm of the adjusted p-value (significance level) for each gene. Genes with high significance and large fold changes are located towards the top and sides of the plot, respectively.
Figure 1. Identification of DEGs data from GEO dataset GSE7116. Volcano plot representing differential gene expression analysis between condition A and condition B. Genes that are significantly up-regulated (log2 fold change >1, adjusted p-value <0.05) are shown in red, while genes that are significantly down-regulated (log2 fold change < −1, adjusted p-value < 0.05) are shown in blue. The x-axis represents the log2 fold change in gene expression between BRONJ vs. control, while the y-axis represents the negative logarithm of the adjusted p-value (significance level) for each gene. Genes with high significance and large fold changes are located towards the top and sides of the plot, respectively.
Ijms 24 08635 g001
Ijms 24 08635 g002 550
Figure 2. GO and KEGG analysis of the DEGs according to (A) biological process, (B) cellular components, (C) molecular functions, and (D) KEGG.
Figure 2. GO and KEGG analysis of the DEGs according to (A) biological process, (B) cellular components, (C) molecular functions, and (D) KEGG.
Ijms 24 08635 g002
Ijms 24 08635 g003 550
Figure 3. Displays a network analysis of the top 150 proteins based on their network parameters: (A) MCC, (B) MNC, (C) Closeness, (D) Betweenness, and (E) Degree, as well as a Venn diagram for hub genes. The nodes in the network represent proteins, and the edges represent their interactions. Panel (F) shows a Venn diagram of the hub genes, with the number of proteins in each intersection indicated. This figure provides insight into the identification of key proteins that may be involved in the network’s function and highlights their importance based on different network parameters.
Figure 3. Displays a network analysis of the top 150 proteins based on their network parameters: (A) MCC, (B) MNC, (C) Closeness, (D) Betweenness, and (E) Degree, as well as a Venn diagram for hub genes. The nodes in the network represent proteins, and the edges represent their interactions. Panel (F) shows a Venn diagram of the hub genes, with the number of proteins in each intersection indicated. This figure provides insight into the identification of key proteins that may be involved in the network’s function and highlights their importance based on different network parameters.
Ijms 24 08635 g003
Ijms 24 08635 g004 550
Figure 4. Displays a network analysis of hub genes, their degree, and co-expression. (A) shows the hub genes in the network. The nodes in the network represent genes/proteins, and the edges represent their interactions. (B) shows a bar chart of the degree of hub genes, which is the number of nodes and edges they interact with. The x-axis shows the number of interactions, and the y-axis shows the hub genes. (C) shows the co-expression of hub genes, which is the level of similarity in expression patterns between genes. This figure provides insight into the key genes that may be involved in the network’s function, their degree of interaction, and their co-expression patterns.
Figure 4. Displays a network analysis of hub genes, their degree, and co-expression. (A) shows the hub genes in the network. The nodes in the network represent genes/proteins, and the edges represent their interactions. (B) shows a bar chart of the degree of hub genes, which is the number of nodes and edges they interact with. The x-axis shows the number of interactions, and the y-axis shows the hub genes. (C) shows the co-expression of hub genes, which is the level of similarity in expression patterns between genes. This figure provides insight into the key genes that may be involved in the network’s function, their degree of interaction, and their co-expression patterns.
Ijms 24 08635 g004
Ijms 24 08635 g005 550
Figure 5. Chemical structures of six compounds: (A) 1,4-dihydronicotinamide adenine dinucleotide, (B) Tanshinone IIA, (C) Glycitein, (D) Napa-bucasin, (E) 3-(5-(4-(Cyclopentyloxy)-2-hydroxybenzoyl)-2-((3-hydroxybenzo[d]isoxazol-6-yl) methoxy) phenyl) propanoic acid, and (F) Alantolactone. The stereochemistry of the compounds is shown using wedges and dashes to indicate the three-dimensional arrangement of atoms in space. Chiral centers are indicated by R or S configurations, and the stereochemistry of specific positions in the molecule is indicated by dashed or solid wedges. The figure also highlights the active site of each molecule, which is an important feature for its function, and this can vary depending on the specific molecule and its interactions with other molecules.
Figure 5. Chemical structures of six compounds: (A) 1,4-dihydronicotinamide adenine dinucleotide, (B) Tanshinone IIA, (C) Glycitein, (D) Napa-bucasin, (E) 3-(5-(4-(Cyclopentyloxy)-2-hydroxybenzoyl)-2-((3-hydroxybenzo[d]isoxazol-6-yl) methoxy) phenyl) propanoic acid, and (F) Alantolactone. The stereochemistry of the compounds is shown using wedges and dashes to indicate the three-dimensional arrangement of atoms in space. Chiral centers are indicated by R or S configurations, and the stereochemistry of specific positions in the molecule is indicated by dashed or solid wedges. The figure also highlights the active site of each molecule, which is an important feature for its function, and this can vary depending on the specific molecule and its interactions with other molecules.
Ijms 24 08635 g005
Ijms 24 08635 g006 550
Figure 6. Heat map displaying the binding affinities (ΔG) of drugs against a panel of genes. The x-axis shows the genes analyzed, while the y-axis displays the drug molecules tested. The color scale represents the affinity range from −4.00 to −11.00 kcal/mol, where blue indicates higher affinity (more negative ΔG values), and red indicates lower affinity (less negative ΔG values). The binding affinities were calculated using molecular docking. (72724: Alantolactone, 439153: 1,4-dihydronicotinamide adenine dinucleotide, 5317750: Glycitein, 10331844: Napabucasin, 23626877: 3-(5-(4-(Cyclopentyloxy)-2-hydroxybenzoyl)-2-((3-hydroxybenzo[d]isoxazol-6-yl) methoxy) phenyl) propanoic acid, 164676: Tanshinone IIA).
Figure 6. Heat map displaying the binding affinities (ΔG) of drugs against a panel of genes. The x-axis shows the genes analyzed, while the y-axis displays the drug molecules tested. The color scale represents the affinity range from −4.00 to −11.00 kcal/mol, where blue indicates higher affinity (more negative ΔG values), and red indicates lower affinity (less negative ΔG values). The binding affinities were calculated using molecular docking. (72724: Alantolactone, 439153: 1,4-dihydronicotinamide adenine dinucleotide, 5317750: Glycitein, 10331844: Napabucasin, 23626877: 3-(5-(4-(Cyclopentyloxy)-2-hydroxybenzoyl)-2-((3-hydroxybenzo[d]isoxazol-6-yl) methoxy) phenyl) propanoic acid, 164676: Tanshinone IIA).
Ijms 24 08635 g006
Ijms 24 08635 g007 550
Figure 7. Docking of small molecule drugs with targets. The active site of protein is shown. Hydrogen bonds and van der Waals interactions between the drugs and protein X are depicted by orange and green dotted lines, respectively. Several important binding interactions are labeled. The molecular docking analysis was performed using AutoDock Vina.
Figure 7. Docking of small molecule drugs with targets. The active site of protein is shown. Hydrogen bonds and van der Waals interactions between the drugs and protein X are depicted by orange and green dotted lines, respectively. Several important binding interactions are labeled. The molecular docking analysis was performed using AutoDock Vina.
Ijms 24 08635 g007
Ijms 24 08635 g008 550
Figure 8. The IL-6 activated by inflammatory triggers initiates the pathway by phosphorylation of STAT3.
Figure 8. The IL-6 activated by inflammatory triggers initiates the pathway by phosphorylation of STAT3.
Ijms 24 08635 g008
Ijms 24 08635 g009 550
Figure 9. General overview of the research pipeline.
Figure 9. General overview of the research pipeline.
Ijms 24 08635 g009
Table
Table 1. Summary of studies that have used GSE7116 dataset.
Table 1. Summary of studies that have used GSE7116 dataset.
StudyStudy AreaStudy ObjectivesNormalizationScreening of Small DrugMolecular DockingReference
Kim et al., 2017BRONJTo identify combined biomarkers associated with BRONJNo preprocessing of dataNONO[13]
Dong Leng et al., 2018Multiple Myeloma (MM)To investigate transcriptional changes of CRISP3 for MM markerCONOR_1.0.1NONO[16]
Jiangnan Sun et al., 2015BRONJTo explore molecular mechanisms associated with BRONJ in patients with MMRMANONO[17]
Juncheng et al., 2020Immune Responses in BRONJTo investigate the immune cellular and genomic profiles of BRONJ and excavate potential small molecule drugNot Mentioned but screened only immune related genesCMapNO[18]
Huanzhi Ma et al., 2021OsteonecrosisTo investigate the DEGs of normal vs. osteonecrosis patients (GSE74089, GSE7116, GSE123568)lmFit and eBayesNONO[19]
Jinpeng Zhuang et al., 2022BRONJTo identify drugs that potentially modulate the risk of BRONJ in cancerGEO2RDGIdbNO[14]
Our Present StudyBRONJTo identify biomarkers and small drug molecules related to BRONJRMACMapYES
Table
Table 2. Details of GEO BRONJ data.
Table 2. Details of GEO BRONJ data.
ReferenceGEO NumberPlatformBRONJ (n)Control (n)RegionReference
Raje et al. (2008)GSE7116GPL5701110USA[40]
Disclaimer/Publisher’s Note: The statements, opinions and data contained in all publications are solely those of the individual author(s) and contributor(s) and not of MDPI and/or the editor(s). MDPI and/or the editor(s) disclaim responsibility for any injury to people or property resulting from any ideas, methods, instructions or products referred to in the content.

Share and Cite

MDPI and ACS Style

Balachandran, K.; Ramli, R.; Karsani, S.A.; Abdul Rahman, M. Identification of Potential Biomarkers and Small Molecule Drugs for Bisphosphonate-Related Osteonecrosis of the Jaw (BRONJ): An Integrated Bioinformatics Study Using Big Data. Int. J. Mol. Sci. 2023, 24, 8635. https://doi.org/10.3390/ijms24108635

AMA Style

Balachandran K, Ramli R, Karsani SA, Abdul Rahman M. Identification of Potential Biomarkers and Small Molecule Drugs for Bisphosphonate-Related Osteonecrosis of the Jaw (BRONJ): An Integrated Bioinformatics Study Using Big Data. International Journal of Molecular Sciences. 2023; 24(10):8635. https://doi.org/10.3390/ijms24108635

Chicago/Turabian Style

Balachandran, Kumarendran, Roszalina Ramli, Saiful Anuar Karsani, and Mariati Abdul Rahman. 2023. "Identification of Potential Biomarkers and Small Molecule Drugs for Bisphosphonate-Related Osteonecrosis of the Jaw (BRONJ): An Integrated Bioinformatics Study Using Big Data" International Journal of Molecular Sciences 24, no. 10: 8635. https://doi.org/10.3390/ijms24108635

APA Style

Balachandran, K., Ramli, R., Karsani, S. A., & Abdul Rahman, M. (2023). Identification of Potential Biomarkers and Small Molecule Drugs for Bisphosphonate-Related Osteonecrosis of the Jaw (BRONJ): An Integrated Bioinformatics Study Using Big Data. International Journal of Molecular Sciences, 24(10), 8635. https://doi.org/10.3390/ijms24108635

Note that from the first issue of 2016, this journal uses article numbers instead of page numbers. See further details here.

Article Metrics

Back to TopTop