0% READ
FULL TEXT
INTRODUCTION
Oral Squamous Cell Carcinoma (OSCC), the most prevalent head and neck tumour, is a major threat to global health. This represents about 90% of OSCC cases, with a projected annual death toll of 177,000 and 350,000 new cases emerging every year (Bray et al., 2018; Sung et al., 2021). The prognosis for OSCC persists to be unfavourable since it is usually aggressive and is often diagnosed in its later stages. Conventional therapies, including radiation, chemotherapy, and surgery, lack efficacy and are typically associated with adverse effects (Aswathy et al., 2024). Although there has been significant progress in the field, its therapeutic outcomes continue to be inadequate. Despite undergoing surgery, 21% - 47% of patients exhibit signs of recurrence (Alim et al., 2024). Patients who had undergone salvage surgery endured a second recurrence within an average of 17 months (Contrera et al., 2022). Despite the fact that adjuvant chemoradiotherapy supplements disease control, it can also cause such adverse effects as trismus, xerostomia, osteoradionecrosis, mucositis, and neck fibrosis (Hosni et al., 2017). These problems underscore the urgent need of low-risk therapeutic intervention, which is patient-oriented.
In this regard, investigators are exploring bioactive phytochemicals derived from medicinal herbs in preclinical trials. For instance, azadirachtin and nimbolide from neem promote apoptosis and downregulate Akt and STAT3 signalling in OSCC cell lines (Prakash et al., 2021). Epigallocatechin-3-Gallate (ECCG), a derivative from Camellia sinensis (green tea), which impeded the proliferation of OSCC cell lines CAL-27 and SCC-25 by inducing S & G₂/M phase cell cycle arrest (Liu et al., 2011). In parallel with such findings in conventional medicine, Homeopathy persists to be a broadly adopted complementary approach of healthcare which has a historical and traditional significance. Its incorporation into modern research-based practice is increasingly being researched through ongoing research aimed at understanding its mechanisms and potential therapeutic utility and efficacy. Galium aparine (G. aparine) is one such medicinal plant that has been used regularly over the years in the field of homeopathy to treat dermal diseases, swollen lymph nodes, and kidney problems (Al-Snafi, 2018; Petkova et al., 2025). Besides the traditional applications, the current studies show that it can have anti-cancer, antioxidant, detoxifying, hepatoprotective, and anti-bacterial effects (Beirami et al., 2024; Bokhari Jasia et al., 2013; Rajkumar et al., 2024; Sahin et al., 2022). Research by Aslanturk et al. reported a study carried out in vitro, which discussed cytotoxicity and induced apoptosis caused by extracts of G. aparine in human breast carcinoma cell lines (MCF 7) and the human Colorectal adenocarcinoma cell line (Caco-2), hence explaining its curative potential (Aslantürk et al., 2017). The G. aparine methanol extract was also tested by another study on breast cancer cell lines, MCF-7 and MDA-MB-231 (Atmaca et al., 2016). It had been found that the extract inhibited MCF-7 and MDA-MB-231 cells without harming non-tumour normal cells from the MCF-10A cell line. These observations indicate that G. aparine could specifically target tumour cells, sparing the healthy cells. Due to its conventional use in treating a variety of ailments, chronicled anti-cancer properties, and corroborating evidence from prior in vitro research, G. aparine has been selected as the herb of preference for subsequent analysis in this study.
Network pharmacology is a promising means to develop medicinal drugs. It employs systems biology and computational techniques to quickly evaluate diseases, targets, and bioactive compounds (Gayathiri et al., 2024; Hopkins, 2008). This integrative technique assists investigators to understand the interactions in biological systems, that guides in identifying vital molecular targets and their signalling pathways associated with the progression of disease. It has been an extensively applied, economic, and high-performing method to study OSCC (Yin et al., 2025; Zhang et al., 2021; Zhao et al., 2024). It also highlights how multi-target bioactive compounds can be used in treating OSCC.
Research objective is to recognize potential bioactive compounds from Gallium aparine that could potentially be utilized for treating OSCC and use computational methods to analyse the multi-target potency of these phytochemicals against the targets identified. In addition, the research seeks to analyse drug-likeness attributes of these compounds through predictive assessments.
MATERIALS AND METHODS
Bioactive Compound Selection
Based on literature evidence, it has been reported that the selected medicinal plant, G. aparine, exhibits potential anticancer action in both in vitro testing and in vivo studies. Indian Medicinal Plant, Phytochemical and Therapeutics (IMPPAT) database https://cb.imsc.res.in/imppat/ has been utilized to retrieve phytochemicals from G. aparine. Identical and duplicate entries were eliminated after the retrieval process. Subsequently, Compounds have been selected according to specified inclusion and exclusion criteria. Inclusion criteria: Alkaloid derivatives, terpenoids, quinazolines, Fatty acids with anticancer potential, iridoids, etc. Exclusion criteria: Sugar, Cholesterol, Alkane, Long chain fatty acids, and Ketone (Vivek-Ananth et al., 2023).
DrugLikeness and Pharmacokinetics
The SMILES (Simplified Molecular Input Line Entry System) representations of the phytochemicals were obtained from the PubChem database (https://pubchem.ncbi.nlm.nih.gov/). These phytochemicals were evaluated for conformity to Lipinski’s rule of 5 to measure drug-like appearance. SwissADME (http://www.swissadme.ch/) web server was used to obtain data for drug-likeness and pharmacokinetics by submitting the SMILES notation of each compound (Daina et al., 2017).
Toxicity Prediction
Lethal dose 50, Carcinogenicity, Acute toxicity classes, Hepatotoxicity, Mutagenicity, and Cytotoxicity were evaluated for a comprehensive assessment of the phytochemicals. The selected ligands were assessed using the ProTox II web server (https://tox.charite.de/protox3/). Phytochemicals with toxic endpoints were excluded (Banerjee et al., 2018).
Biological Activity Spectrum (BAS) Analysis
To predict the biological activities and pharmacological effects of the lead compounds, the Prediction of Activity Spectra for Substances (PASS) (http://www.way2drug.com/passonline/predict.php) online web server has been used. Server operates through input of the canonical SMILES notations of each phytochemical and detects the probability of the compound to possess certain biological effect. These predictions enable the understanding of the pharmacological properties of phytochemicals.
OSCC Target Analysis
Genes associated with OSCC were gathered from three databases: Geo, DisGeNet, and GeneCards, which facilitated the aggregation and compilation of data on OSCC-related genes. The keywords employed for this investigation encompassed "Oral cancer." “OSCC" and "Head and Neck Cancer." A Venn diagram was created using Interactive Venn (https://www.interactivenn.net/) to demonstrate intersection of these three data sets.
A compound-disease-target analysis was conducted to discover possible targets common to the phytochemicals and OSCC. This method finds the targets of bioactive chemicals and juxtaposes them with recognized targets associated with OSCC.
PPI Network Construction and Hub Gene Identification
STRING database (https://string-db.org/) was utilized for further analysis of genes obtained from compound-disease-target-analysis. The species of interest was given “Homo sapiens,” and a confidence score of 0.400 (medium confidence) was considered. PPI network has been generated and visualized by using Cytoscape version 3.10.0 where hub gene identification was carried out via the CytoHubba plugin. The hub genes were detected using topological analysis methods, including degree and closeness centrality based on node ranking (Szklarczyk et al., 2023).
Gene ontology and KEGG
Gene Ontology (GO) and Kyoto Encyclopaedia of Genes and Genomes (KEGG) evaluation has been shown to elucidate functional roles and pathway involvement of target genes MMP1 and MMP9 in OSCC. GO enrichment—covering Molecular Functions (MF), Biological Processes (BP), and Cellular Components (CC) — was performed utilizing the Enrichr web server. Official gene symbols for (Matrix Metalloproteinase-1) MMP1 and (Matrix Metalloproteinase-9) MMP9 were submitted, and enriched GO terms were reviewed via bar plots generated by Enrichr. KEGG pathway enrichment was also performed using Enrichr (https://maayanlab.cloud/Enrichr/), with “Homo sapiens” selected as the species. Enriched pathways, including “Pathways in Cancer,” were identified and visualized through bar plots for further interpretation (Chen et al., 2013; Kuleshov et al., 2016; Xie et al., 2021).
Expression, Correlation, and Survival Profile Analysis of Hub Genes
Gene Profiling Interactive Analysis (http://gepia.cancer-pku.cn/) (GEPIA), having a beneficial network tool analysing expression and correlation data of the hub genes. The tumour data from the Genotype-The Cancer Genome Atlas (TCGA) and Genotype Tissue Expression Project (GTEx) have been compiled and employed by GEPIA. Survival analysis was carried out by means of GEPIA to analyse the prognostic relevance of gene expression patterns in OSCC (Tang et al., 2017).
Correlation between Target Genes and Infiltration of Immune Cells
To evaluate the immune cell infiltration within tumours, the web-based TIMER database (http://timer.cistrome.org/) was employed with default settings. Correlation between target genes and infiltration of immune cells. Tumour Immune Estimation Resource (TIMER) was used to investigate the interrelation between target genes and immune cells in tumour tissues. Subsequently, the correlation between the expression levels of the two hub genes and immune cell infiltration was studied using the TIMER database. Correlation of immune cell inflammation (including T-Cell CD8+, B cell, and T-Cell CD4+) and hub genes in patients with Head and Neck Squamous Cell carcinoma (HNSC) was performed. Spearman's rank correlation was utilized in Timer to determine the relationship between gene and the level of immune infiltration. A p-value < 0.05 was considered statistically significant (Li et al., 2017).
Tissue presence/cell line analysis
We have evaluated the expression of the selected gene targets in several human tissues and cell lines using the findings of the ProteomicsDB platform (https://www.proteomicsdb.org/), This protein expression analysis enabled us to explain the potential activity of the compound in some tissues and cell lines, revealing some possibilities in terms of therapeutic implications (Samaras et al., 2020; Schmidt et al., 2018).
Hub Gene Druggability Analysis
The Druggability of the hub genes has been assessed in the Drug Gene Interaction Database (DGIdb) (https://www.dgidb.org/). It is a positive site that unites information on drug-gene interactions based on several repositories and literature sources. The tool gives data regarding known interactions between medications and genes. The strategy enables a comprehensive analysis of the druggability of hub genes and helps to identify new therapeutic interventions of the hub genes (Cannon et al., 2024).
Molecular docking
Docking analysis was conducted to determine the interactions of the bioactive phytochemicals with the identified targets and thus confirm network pharmacology predictions. The approach provides a considerable amount of information regarding the multitargeting action of phytochemicals that simultaneously act on different proteins. This enables the realization of the possible interactions and synergistic effects between compounds and targets.
RCSB PDB (https://www.rcsb.org/) has been used to download 3D structures of target proteins, and the PubChem database was used to download 3D structures of phytochemicals. The PDB ID for MMP1 is 3SHI, and MMP9 is 1ITV. BIOVIA Discovery Studio Visualizer v21.1.0.20298 has been employed for preparing the target protein. At first, to eliminate solvent interference, all water molecules adjacent to the protein were removed. The ligand was deleted to expose the active site, after the co-crystallized ligand coordinates were recorded. Nonpolar hydrogens were merged, Polar hydrogens were included, and structure of the protein was held flexible during docking (Berman et al., 2000).
To evaluate the phytochemicals' binding affinities, PyRx was used for molecular docking. Compounds were treated as flexible particles by enabling all rotatable bonds to move freely. This facilitated the ability to identify the strongest binding conformations. A grid box was defined with dimensions of 11.52 Å (x-axis), -28.09 Å (y-axis), and 6.77Å (z-axis) for MMP1 and -42.21 Å (x-axis), -30.62 Å (y-axis), and -7.27 Å (z-axis) for MMP9, focused on the critical amino acid residues essential for ligand interaction. Out of the ten predicted binding configurations, the conformation displaying the lowest value was chosen as being the most likely pose. This protein-ligand complex was evaluated employing BIOVIA Discovery Studio Visualizer v21.1.0.20298.
Swiss Similarity
Phytochemicals with significant binding affinities, favourable drug-likeness properties, and that met the required criteria for toxicity analysis were subsequently analysed utilizing the Swiss similarity (http://www.swisssimilarity.ch/) web server to recognize compounds with prior FDA authorization (Bragina et al., 2022).
Bioavailability Radar of the leads
A webserver tool, SwissADME (http://www.swissadme.ch/), served to assess the bioavailability radar of the phytochemicals, promptly demonstrating if the compounds have oral bioavailability.
Mutation Analysis
The mutation rates of the target genes MMP1 and MMP9 in HNSC were analyzed using the cBioPortal (https://www.cbioportal.org) platform. The dataset selected for this analysis was Head and Neck Squamous Cell Carcinoma (TCGA, PanCancer Atlas), comprising 523 samples. Both targets were evaluated for the frequency and types of genetic alterations present within this cohort (de Bruijn et al., 2023).
In silico Mutagenesis
Since the 3D structure of the mutated form of MMP1 is not available in the Protein Data Bank (PDB), the mutation was introduced computationally. Among them, Asn180—a residue of the active site—plays a supportive role in maintaining the binding site configuration. To investigate its influence on ligand interaction, Asn180 was mutated to alanine (N180A) using PyMOL software. The 3D structure of the mutated form of MMP9 is available in the Protein Data Bank under the PDB ID: 2OVX.
Physiochemical Characterization
The mutated MMP1 structure was validated using the SAVES server (https://saves.mbi.ucla.edu/), utilizing ERRAT, Verify 3D, and PROCHECK tools (Colovos & Yeates, 1993; Pontius et al., 1996).
Analysis of dynamic simulations
The GROMACS 2022.3 software suite with the Amber force field was used to conduct Molecular Dynamics (MD) simulations using all protein topologies. Crystal structures of MMP1wild and MMP1mutant, MMP9wild and MMP9mutant were retrieved in the Protein Data Bank, and the protein-ligand complexes were constructed in the presence of solvent and ions to simulate the physiological conditions. Delta-7 Avenesterol ligand was parameterized, and its topology was generated through the ACPYPE toolkit with GAFF compatibility, which ensures that it is properly integrated with the Amber force field environment. The solvation of each system was done in a periodic cubic water box using TIP3P water models, and counter-ions were added where needed. The steepest descent algorithm was used to minimize energy and remove steric clashes, and optimize the initial geometry. After the minimization, the first equilibration was obtained in two phases: Bussi-velocity rescaling thermostat with a 100 ps NVT (constant Number, Volume, Temperature) ensemble to stabilize the temperature and Parrinello-Rahman barostat with a 100 ps NPT (constant Number, Pressure, Temperature) ensemble to stabilize the system pressure and density. Position restraints were applied to heavy atoms during equilibration. The production MD simulation was then run for 100 ns under the NPT ensemble with periodic boundary conditions. A time step of 2 fs was used, constraining all bonds involving hydrogens by the LINCS algorithm. Long-range electrostatics were treated with the Particle Mesh Ewald (PME) method, and van der Waals interactions were truncated at 1.0 nm. Trajectory analyses, including calculations of RMSD, RMSF, protein-ligand RMSD, ligand RMSD, radius of gyration, and hydrogen bonding, were performed using standard GROMACS utilities. Visualizations and binding pocket analyses were conducted using PyMOL (Berendsen et al., 1995).
RT-PCR Analysis of MMP1 and MMP9
KB oral carcinoma cells were treated with the IC₅₀ concentration ofG. aparine extract (272.1 µg/mL), determined from preliminary cytotoxicity screening, and incubated for 24 h. Untreated cells maintained under identical culture conditions served as controls. Following treatment, cells were harvested for total RNA isolation.
The standard procedures involved the extraction of total RNA with TRIzol ™ reagent. The concentration and purity of the RNA were determined using a NanoDrop spectrophotometer, and samples that had acceptable A260/A280 ratios were subjected to downstream analysis. A 1 µg of total RNA was synthesized into complementary DNA (cDNA) using Bio-Rad iScript™ cDNA Synthesis Kit according to the instructions of the manufacturer.
iTaq™ Universal SYBR® Green Supermix was used on a Bio-Rad CFX Opus real-time PCR system to do quantitative real-time PCR. MMP1, MMP9 were subjected to gene-specific primers, and β-actin was the internal control. All the reactions were repeated. The 2⁻ΔΔCt method was used to obtain relative expression as a fold change compared to untreated control cells.
Primer Sequences
RESULTS
Bioactive Compounds
The IMPPAT database retrieved 106 bioactive compounds that were reported in G. aparine plants. The inclusion and exclusion criteria were used to select 38 phytochemicals to be further analyzed. The list of compounds with PubChem ID and 2D structures is represented in the form of Supplementary Table S1.
| Forward | Reverse | |
|---|---|---|
| β-actin | AGA GCT ACG AGC TGC CTG AC | AGC ACT GTG TTG GCG TAC AG |
| MMP1 | AGC CTG GTG ATG AGG ACT GA | TGG GAT GGT CTT CTC CAG GA |
| MMP9 | CCT GGA GAC CTG AGA ACC AA TC | GAT TTC GAC TCT CCA CGC AT CT |
DrugLikeness
The 38 phytochemicals of G. aparine were tested on their drug-likeness. One of them, Asperulosidic acid and Deacetylasperulosidic acid, had two Lipinski violations, and another, Rutin, had three violations; therefore, these compounds were excluded from further pharmacokinetic study. Asperulosidic acid contains 6 hydrogen bond donors, 12 acceptors, a low Log P of -2.02, and a high TPSA of 192.44 Ų. Rutin has 16 hydrogen bond acceptors, 10 hydrogen bond donors, a Log P of -1.51, and a molecular weight of over 500 Da. Such high parameters can reduce bioavailability through the mouth. Conversely, 35 compounds with one violation each were selected for further research. In general, most of the phytochemicals met the necessary drug-likeness properties, which implies their suitability in future pharmacokinetic studies. Supplementary Table S2 reflects the drug-likeness of the phytochemicals.
Pharmacokinetics and Toxicity
The pharmacokinetic analysis of the 35 phytochemicals used showed that 19 of them were Blood-Brain Barrier (BBB) permeant and thus not subjected to further analysis. Among the other 16 compounds that were not permeable to the BBB, the majority of them displayed low Gastrointestinal Absorption (GIA), and this may potentially restrict oral bioavailability. Neophytadiene was found to be one of these substrates of P-glycoprotein (P-gp). This resulted in the exclusion of 20 compounds in the further toxicity screening, and the remaining 15 compounds were to be under further study. Supplementary Table S3 is a table that shows the pharmacokinetics of the phytochemicals.
On the toxicity screening of the retained subset, such factors as LD50 values, toxicity classes, carcinogenicity, hepatotoxicity, mutagenicity, and cytotoxicity were taken into consideration.
Toxicity screening of the 15 shortlisted compounds demonstrated that Luteolin was likely to be carcinogenic, whereas Calamenene was likely to be mutagenic; neither was to be further assessed. The rest 13 compounds were categorized into 4 to 6 toxicity classes, which are moderate and low toxicity. Their LD₅₀ values ranged from 890 mg/kg (Campesterol, Stigmasterol, and β-Sitosterol) to 70,000 mg/kg (β-Amyrin). Importantly, none of these 13 compounds showed any significant predicted toxicity, supporting their suitability for further investigation. Supplementary Table S4 represents the toxicity analysis of the 15 phytochemicals.
Biological Activity Spectrum (BAS) Analysis
To predict the biological actions of 13 compounds using the PASS server, the analysis revealed that every compound exhibited antineoplastic activity (Supplementary Table S5). The p-value was between 0.249 and 0,916. The fact that the Probability of activity (Pa) exceeds the Probability of inactivity (Pi) implies that they might also produce the anticipated antineoplastic effect.
OSCC Target Identification
Data on genes from 3 public databases were mined to extract genes associated with OSCC. These were the Gene Expression Omnibus (GEO, 143 genes), GeneCards (5050 genes), and DisGeNET (734 genes), which produced a total of 5927 genes.
GEO was accessed to retrieve the GSE37991 dataset, which has 40 samples of the oral squamous cell carcinoma and 40 samples of normal tissues. DEGs were identified with GEO2R with the use of |human|>DEGs were determined using GEO2R using |human| Log2FC| > 2 and adjusted p < 0.05 to identify upregulated and downregulated genes, respectively. A volcano plot was generated to visualize the results (Figure 1a).
From GEO2R, a total of 477 Differentially Expressed Genes (DEGs) were identified, comprising 130 upregulated genes and 347 downregulated genes.
A Venn diagram was developed to visually demonstrate the common genes among GEO, GeneCards, and DisGeNET (Figure 1b). Venn analysis revealed that out of 5,927 genes, a key set of 17 overlapping genes was identified, showing a strong association with OSCC. The 17 shared OSCC associated genes are: EPHB2, MMP1, MPP9, CSF2, SERPINE1, SPP1, CDKN2A, PLAU, CA9, MMP13, PTHLH, CXCL9, PDPN, MMP10, CXCL13, PAEP and FBLIM1.
PPI Network and Hub Gene Identification
A Protein-Protein Interaction (PPI) network was constructed for the 17 overlapping genes using the STRING database with a confidence score threshold of 0.400, revealing significant interactions among the target proteins PPI enhancement p-value of <1.0e-16 indicates that the observed interactions are statistically significant.
Cytoscape (version 3.10.0) has been utilized to visualize the network, revealing 17 nodes and 47 edges. The significant genes in the network were determined by employing degree centrality analysis. MMP1 and MMP9 were recognised as the highly degree genes based on degree and closeness centrality (Figure 1c and 1d). Both these matrix metalloproteinases have been associated with tumour invasion and cancer metastasis, in addition to the degradation of the extracellular matrix. Their central alignment in the PPI network could suggest their involvement as critical mediators of anticancer effects of the analysed compounds.
Gene Ontology And KEGG
Gene ontology analysis of MMP1 and MMP9 demonstrated significant enrichment across key biological processes, molecular functions, and cellular elements.
Biological procedures included pathways relating to extracellular matrix remodelling and cellular stress responses (Figure 2a).
Regarding cellular components, these genes were predominantly associated with extracellular matrix related and granule associated structures (Figure 2b).
In terms of molecular functions, MMP1 or MMP9 showed strong enrichment in proteolytic and metal ion binding activities consistent with their enzymatic roles (Figure 2c).
Our examination of principal pathways connected to two target genes, MMP1 and MMP9, in HNSC revealed notable pathway enrichment. Specifically, both genes were significantly enriched in the “Pathways in Cancer” KEGG category, with a p-value of 7.04E-04 indicating strong statistical significance (Figure 2d).
Expression, Correlation, and Survival Profile Analysis of Hub Genes
Figure 3a demonstrates box plots of hub genes’ expression profiles in OSCC compared with normal tissue samples. It is revealed from the box plot that MMP1 and MMP9 have exhibited increased expression levels in OSCC samples than in normal tissue samples. The height of the bar depicts the median expression level of normal tissue and the tumour sample. The whiskers arising from the bar and extending out relate to the range of expression, as well as the maximum and minimum values for every sample set.
The correlation between the hub genes in OSCC has been examined using a correlation study. The results revealed that MMP1 and MMP9 were positively correlated, suggesting a moderate association between the expression of these target genes. This provides insights into potential interactions and mechanisms of these genes, demonstrating their significance in OSCC development and progression. The network of the hub genes displaying positive correlation in OSCC is illustrated in Figure 3b.
Data retrieved from GEPIA revealed that elevated MMP1 expression levels were connected to worse overall survival in HNSC (HR = 1.3, P = 0.047). On the contrary, the expression levels of MMP9 were statistically insignificant. (HR 0.99, P = 0.9) Figure 3c.
Correlation between Target Genes and Infiltration of Immune Cells
The correlation analysis results revealed a notable association between the expression of MMP1 and infiltration of immune cells in HNSC. Specifically, MMP1 expression showed positive correlation with the infiltration levels of CD8+ T cells, B cells, and CD4+ T cells, suggesting its potential role in modulating the tumour immune microenvironment. Conversely, MMP9 had a different pattern: although its expression had a positive correlation with the infiltration of CD4+ T cells, it had a negative correlation with the infiltration of CD8+ T cells and B cells. Such results suggest that both MMP1 and MMP9 might have a role in the recruitment of immune cells and the HNSC tumour microenvironment, and may have implications on immune evasion and tumour progression (Figure 4a & b).
Tissue presence/cell line analysis
MMP1 and MMP9 display distinct expression patterns in various human tissues as well as cell lines. In order to more effectively comprehend and analyze these expression profiles, a heatmap was generated from ProteomicsDB (Figure 5), which pictorially illustrates the relative concentration of both genes using a color gradient scale. Each shade in the scale corresponds to a specific degree of expression, providing a straightforward approach to analyzing the data. Findings reveal that, scope of Neck and head squamous carcinoma, MMP1 displayed significant protein expression within oral epithelium, saliva, and salivary gland tissues; likewise, MMP9 exhibited marked expression in the BHY cell line (OSCC cell line) and PCI 6A (HNSC).
Druggability of targets
Druggability evaluation showed target genes linked to drug interactions, particularly focused on the interaction score threshold of 0.17 and above. Among the studied genes, MMP9 presented the maximum amount of drug interactions within the set limit compared to MMP1. MMP1 interacted with 20 drugs, while MMP9 interacted with 25 drugs (Table 1).
| Sl. No. | Target | Gene |
|---|---|---|
| 1 | MMP1 | Leflunomide, therapeutic glucocorticoid, collagenase clostridium histolyticum marimastat, rebimastat, hydrocortisone butyrate, doxycycline calcium, medroxyprogesterone, interferon beta, doxycycline anhydrous, cobalt, vitamin E, pentosan polysulfate sodium, lamivudine, sirolimus, leuprolide acetate, tyrosine kinase inhibitor, ribavirin, triamcinolone, immunosuppressant. |
| 2 | MMP9 | Curcumin pyrazole, tst10088, celecoxib, curcumin, tozuleristide, shark cartilage extract ae-941, hydralazine, ulinastatin, bevacizumab, carboxylated glucosamine, methyldopa anhydrous, and ecaliximab, dp-b99, neovastat, prinomastat, rebimastat, incyclinide, cts-1027, demethylwedelolactone, s-3304, ilomastat, marimastat, azd-1236, antihypertensive agent, nifedipine. |
Molecular docking
Docking was carried out between the two targets, MMP1 and MMP9, with 38 phytochemicals. Phytochemicals with binding affinity falling within the range of -6.5 to -8.0 kcal/mol were considered to have stronger interactions with the target proteins (Figure 6a and 6b).
The binding affinities of selected phytochemicals against target MPP1 ranged between -6.5 to -7.8 Kcal/mol. The binding affinities of all 13 phytocompounds fell within the range of -6.5 to -8.0kcal/mol for MMP1. Notably, the control drug doxorubicin exhibited a binding affinity of -7.3 kcal/mol toward MMP1, indicating that several of the tested compounds demonstrated comparable or stronger binding potential than the control drug. Table 2 represents the binding affinities of phytochemicals with MMP1.
| Sl. No. | Compound | Binding affinity | No Of H Bond | Conventional hydrogen bond | Carbon Hydrogen / Pi-Anion /Pi-Sulfur /Pi-Pi Stacked / Pi- Alkyl Bond | Van Der Waals |
|---|---|---|---|---|---|---|
| 1 | Beta-Amyrin | -6.7 | 0 | - | His218, Leu181 | His222, Ala184, Glu219, His183, Ala182, Val215, Tyr240, Ser239, Gly179, Pro238, Asn180, His228 |
| 2 | Chlorogenic Acid | -7.5 | 5 | Pro173, Glu219, Ala182, Gly179, Leu181 | Pro238, Asn180 | Ser172, Phe174, His168, His183, His218, Val215, Tyr240, Ser239 |
| 3 | Cycloartenol | -7 | 0 | - | Phe185, His218, Pro238 | Gly179, Ser239, Leu181, Ala182, Asn180, Tyr240, His183, Glu219, His222, His228, Ala184, Ser172 |
| 4 | Asperuloside | -7.8 | 4 | Glu219, Ala184, Asn180, Arg214 | His183 | Leu235, Val215, Tyr237, Thr241, Tyr240, His218, Ser239, Leu181, Ala182, Gly179, Pro238, His228, His222 |
| 5 | Alpha- Curcumene | -6.5 | 0 | - | Val215, Glu219, Tyr240, Ser239, His218, Leu181, Tyr210 | Arg214, Pro238, Asn180, Gly179, Ala182 |
| 6 | 15-Demethyl Plumieride | -6.9 | 5 | Asn180, Gly179, Leu181, Ala182, His183 | His228 | Phe174, Ser172, Glu219, His218, Val215, His222, Pro238 |
| 7 | Delta7-Avenasterol | -7.2 | 0 | - | His183, Phe185, His228, His222, Leu181 | Glu219, His218, Ala182, Ala184, Tyr240, Ser172, Asn180, Pro238, Ser239, Gly179, Gln186 |
| 8 | Cycloartanol | -7.3 | 0 | - | Leu181, Val215, His218, Tyr240 | Glu219, His222, Gly179, Ala182, Asn180, His228, Ala184, His183, Pro238, Ser172, Arg214, Ser239 |
| 9 | Lanosterol | -6.5 | 0 | - | His218, Val215, Leu181 | Gly179, Ser172, Asn180, His183, Ala184, His222, His228, Ala182, Glu219, Pro238, Ser239, Tyr240 |
| 10 | (3S)-3-Hydroxy-2,3-Dihydro-1H-Pyrrolo [2,1-B]Quinazolin-9-One | -7.2 | 0 | - | Gln154, Glu60, Ala13, Pro62, Arg106, Glu14 | Val58, Phe59, Phe153, Ala104, Met112, Leu105 |
| 11 | Campesterol | -7.3 | 0 | - | Leu181, His218, His222 | Gly179, Glu219, Ala182, His228, Asn180, His183, Ala184, Val215, Tyr240, Ser239, Pro238, Ser172, Phe185 |
| 12 | Beta-Sitosterol | -6.8 | 0 | - | Ser239, Phe185 | Gly179, Ser172, Leu181, Asn180, His228, Tyr240, Pro238, Glu219, Ala182, His218, His222, Ala184, His183 |
| 13 | Stigmasterol | -7.2 | 0 | - | His218, His222, Leu181, Phe185 | His218, Glu219, Ala182, Pro238, Tyr240, His183, Ala184, Asn180, Ser239, Gly179, Gln186, Ser172 |
| C | Doxorubicin | -7.3 | 3 | Glu219, Ala184, Asn180 | Pro238, Ser239 | His218, His222, Leu181, His183, Ala182, Tyr237, Tyr240, Val215, Gly179 |
For MMP9, the binding affinities of the phytochemicals ranged from -4.8 to -7.0 kcal/mol. Among the 13 compounds, (3S)-3-hydroxy-2,3-dihydro-1H-pyrrolo[2,1-b]quinazolin-9-one exhibited the strongest interaction with a binding affinity of -7.0 kcal/mol. In contrast, the reference drug doxorubicin demonstrated a higher binding affinity of -8.2 kcal/mol, indicating stronger binding to the target. Notably, only four phytochemicals—(3S)-3-hydroxy-2,3-dihydro-1H-pyrrolo[2,1-b]quinazolin-9-one, Δ7-Avenasterol, Stigmasterol, and Asperuloside—showed binding affinities within the favorable range, though all remained slightly weaker compared to the control drug. Table 3 represents the binding affinities of phytochemicals with MMP9.
| Sl. No. | Compound | Binding affinity | No of H Bond | Conventional hydrogen bond | Carbon Hydrogen / Pi-Anion /Pi-Sulfur /Pi-Pi Stacked / Pi- Alkyl Bond | Van Der Waals |
|---|---|---|---|---|---|---|
| 1 | Beta-Amyrin | -4.8 | 1 | His150 | Lys54 | Asp164, Ile8, Asp56, Ser71, Gly72, Gln100 |
| 2 | Chlorogenic Acid | -6.4 | 2 | Asp10, Leu55 | Asp164 | Gln163, Asp148, Phe9, His150, Ile8, Lys54, Gln100, Thr102, Asp56 |
| 3 | Cycloartenol | -5.8 | 0 | Thr102, His150, Asp10, Lys23, Asp56, Gln100, Lys54, Ile8, Gly72 | ||
| 4 | Asperuloside | -6.5 | 2 | Gln154, Glu157 | Arg106 | Val170, Phe142, Asn177, Leu113, Leu105, Ala159, Ser172, Ser107, Glu175, Lys158 |
| 5 | alpha-Curcumene | -5 | 0 | Ala159, Leu105, Leu113 | Arg106, Ser107, Gln154, Asn177, Phe142, Ser172, Val170 | |
| 6 | 15-Demethyl plumieride | -6 | 2 | Asp10, Lys54 | Gly72, Phe9 | Ile8, His150, Asp151, Asp56, Arg73, Gln100, Thr102, Leu55, Ser71, Gly117 |
| 7 | Delta7-Avenasterol | -6.8 | 0 | - | Ile8 | Lys23, Asp10, Lys54, Asp56, Gln100, Thr102, His150, Gly117, Asp148, Thr149 |
| 8 | Cycloartanol | -5.8 | 0 | Ile8, Lys23, Lys54, Asp10, His150, Asp56, Gln100, Thr102, Asp148, Thr149, Gly117 | ||
| 9 | Lanosterol | -5.6 | 0 | Ile8 | Lys23, Asp24, Lys54, His150, Asp10, Thr102, Gln100, Asp56 | |
| 10 | (3S)-3-hydroxy-2,3-dihydro-1H-pyrrolo [2,1-b]quinazolin-9-one | -7 | 0 | - | Gln154, Glu60, Ala13, Pro62, Arg106, Glu14 | Val58, Phe59, Phe153, Ala104, Met112, Leu105 |
| 11 | Campesterol | -6.4 | 0 | Ile8 | Asp10, Lys23, Asp56, Gln100, Thr102, His150, Lys54, Gly117, Asp148, Thr149 | |
| 12 | beta-Sitosterol | -6.4 | 0 | Ile8 | Lys23, Asp10, Asp56, His150, Gln100, Asp148, Thr149, Gly117, Thr102, Lys54 | |
| 13 | Stigmasterol | -6.8 | 0 | - | Ile8 | Leu55, Lys23, Asp10, Asp56, His150, Thr102, Thr149, Gln100, Gly117, Asp148, Lys54 |
| C | Doxorubicin | 6 | His150, Asp148, Arg118, Ser116, Gln100, Asp56 | Ala99, Asp97 | Gly163, Asp164, Thr149, Arg119, Val98, Gly117, Thr102, Gly72 |
Swiss Similarity
The identified 4 leads (Asperuloside, Delta7-Avenasterol, (3S)-3-hydroxy-2,3-dihydro-1H-pyrrolo[2,1-b]quinazolin-9-one, and Stigmasterol) with high binding affinity, which were common between both MMP1 and MMP9, were assessed for their structural similarity relative to FDA-authorized drugs by means of SWISS similarity.
The data consisting of FDA-authorized drugs revealing similarity to the selected compounds are provided in Supplementary Table S6 (a-d). Beta Sitosterol is one such compound that exhibited a high similarity score for the compounds, delta7-Avenasterol, Campesterol, and Stigmasterol with similarity scores of 0.655 and 0.850, respectively (Zoete et al., 2016). These observed computational interactions align with results from a study previously published, which revealed that Beta Sitosterol lowered oxidative stress and restored the antioxidant status in DMH-induced colon carcinoma in Wistar rats, thereby serving as a prospective antineoplastic compound for cancer therapy (Baskar et al., 2012). Similarly, Trehalose is one such compound that showed similarity to Asperuloside. In vivo, trehalose liposomes were shown to suppress lung carcinoma growth in xenograft-bearing mice by inducing apoptosis (Ichihara et al., 2017). Likewise, trehalose inhibited malignant melanoma by inducing G2/M cell-cycle arrest and apoptosis, suggesting its potential as a topical therapeutic agent for ocular surface malignancies (Kudo et al., 2012). More recently, trehalose was reported to reduce breast cancer cell growth by modulating nucleotide and NAD⁺ metabolism, with more pronounced effects in cancer cells than in endothelial cells (Opielka et al., 2020). Thiodigalactoside is another compound that showed similarity to Asperuloside. In vivo studies demonstrated that TDG inhibited melanoma and breast cancer growth by blocking galectin-1-mediated tumour-promoting pathways, including angiogenesis, immune dysregulation, and oxidative stress protection (Ito et al., 2011). Similarly, TDG suppressed oral squamous cell carcinoma by targeting β-galactoside-binding protein and inhibiting regulatory T-cell activity, resulting in apoptosis, cell-cycle arrest, and anti-angiogenesis (Aggarwal & Das, 2016).
Bioavailability Radar
The bioavailability radar allows for to quick determination of the drug-likeness of a phytochemical. The plot should lie within the pink space for the phytochemical to be considered to have drug-likeness. Polarity and flexibility represent two key factors that influence the compound’s bioavailability. Assessing polarity through TPSA suggests, phytochemicals having a value greater than 20 Å2 and less than 140 Å2 have higher oral bioavailability. The degree of flexibility depends on rotatable bonds. Phytochemicals possessing more than ten rotatable bonds are likely to show reduced oral bioavailability. Three of the compounds that are investigated meet the radar plot standards; therefore, they are predicted to have good oral bioavailability (Figure 7).
Mutation Analysis
Based on data retrieved from cBioPortal, both genes show mutations in the opted data sources containing 523 samples. MMP1 showed 6% mutation, and MPP9 showed less than 1% mutation (Figure 8a). Co-expression analysis between the target genes MMP1 and MMP9 showed a positive correlation (Figure 8b). Frequency and location of genetic alterations in MMP1 and MMP9 are depicted in Figure 8c. In MMP1, five missense mutations have been reported at amino acid positions 214, 224, 308, 412, and 451, and one nonsense mutation has been reported at amino acid position K356.
In silico Mutagenesis and Docking
The structural validation of the mutated MMP1 model confirmed good overall quality with acceptable stereochemical features. The ERRAT analysis yielded an overall quality factor of 67.347. Ramachandran plot assessment using PROCHECK revealed that 91.6% of residues were located in the most favored regions, while 8.4% fell within the additionally allowed regions. These results indicate that the model possesses a well-defined backbone geometry and is reliable for subsequent structural and docking analyses (Supplementary Figures S1a & b). The Verify3D analysis showed that 98.75% of the residues achieved an average 3D-1D score ≥ 0.1, confirming the reliability of the model structure.
A total of 13 phytocompounds were docked against the active site mutant targets, MMP1 and MMP9, to assess their binding potential.
For mutated MMP1, the binding affinities of the 13 phytochemicals ranged from -6.4 to -8.1 kcal/mol. The reference drug doxorubicin demonstrated a stronger binding affinity of -8.7 kcal/mol. Among the phytochemicals, Campesterol exhibited the highest affinity (-8.1 kcal/mol), while alpha-Curcumene showed the lowest (-6.4 kcal/mol). Table 4 summarizes binding affinity and amino acid interactions with mutated MMP1. Figure 9a shows the 2D interaction of the compounds against the mutated MMP1, respectively.
| Sl. No. | Compounds | Binding Affinity (Kcal/mol) | No Of H Bond | Hydrogen Bond | Carbon Hydrogen / Pi-Anion /Pi-Sulfur /Pi-Pi Stacked / Pi- Alkyl Bond | Van Der Waals |
|---|---|---|---|---|---|---|
| 1 | beta-Amyrin | -7.3 | 1 | Ala182 | Phe185 | Gly192, Gln186, His228, His222, Leu181, Ala180, His183, Ala184, Phe174, Ser172 |
| 2 | Chlorogenic acid | -7.8 | 4 | Ala184, Arg214, Tyr237, Ser239 | His183, Val215, His218, Glu219, Leu181 | Tyr240, Leu235, Thr241, Pro238, His228, His222, Ala182 |
| 3 | Cycloartenol | -7.5 | - | - | His218, Val215 | Ser239, Tyr240, Pro238, Leu181, Ala182, Ala180, His183, Phe174, Ser172, Phe185, His228, Glu219, Ala184, Gln186 |
| 4 | Asperuloside | -7.7 | 2 | Arg214, Ala184 | His183 | Thr241, Leu235, Tyr237, Val215, Tyr240. Ser239, His218, Ala182, Pro238, Gly179, Leu181, His228, Glu219, Ala180, His222 |
| 5 | alpha-Curcumene | -6.4 | - | - | Leu181, Val215, Glu219, Tyr240, Ser239, His218, Tyr210 | Arg214, Ala182, Ala180, Gly179, Pro238 |
| 6 | 15-Demethyl plumieride | -7.1 | 1 | Tyr237 | His218, Pro238, Glu219, Leu181 | Leu235, Arg214, Ser239, Gly179, Tyr240, Ala180, Ala182, Tyr210, Glu209 |
| 7 | delta7-Avenasterol | -7.2 | - | - | Tyr240, Val215, His218, Leu181, Ala180 | His228, His222, Glu219, Pro238, Ser239, His183, Ala182, Phe174, Ser172, Pro173 |
| 8 | Cycloartanol | -7.2 | - | - | Leu181, His218 | Ser239, Tyr240, Ala180, His183, Pro238, Gly179, Ala182, His228, Val215, Glu219, His222, ALa184 |
| 9 | Lanosterol | -7.5 | - | - | Leu181, Val215, His218, Tyr240, Pro173 | Glu219, Pro238, Arg214, Ala182, His183, Thr241, Ala180, Phe174, Ser172, Asp175 |
| 10 | (3S)-3-hydroxy-2,3-dihydro-1H-pyrrolo [2,1-b]quinazolin-9-one | -7.1 | 3 | Pro238, Ala182, Leu181 | Ser239, Val215, Glu219, | Ala180, His228, His218, Tyr240, Arg214 |
| 11 | Campesterol | -8.1 | - | - | His218, Leu181 | Ser172, Ala182, Gly179, Ala180, Ser239, Thr241, Tyr237, Arg214, Tyr240, Val215, Pro238, Glu219, His228, His183, Ala184, His222 |
| 12 | beta-Sitosterol | -7.4 | - | - | His218, Leu181 | Gly179, Thr241, Ser239, Tyr237, Tyr240, Pro238, His228, Ala180, His222, Arg214, Val215, Ala182, His183, Glu219, Phe185, Ser172, Ala184 |
| 13 | Stigmasterol | -7.2 | - | - | Phe185, His228, His222 | Ser239, Leu181, Ala180, Gly179, Tyr240, Pro238, Ala182, His218, Glu219, His183, Ala184, Gln186, Ser172 |
| C | Doxorubicin | -8.7 | 4 | Ser239, Tyr237, His222, Pro238 | Leu181, Tyr219, Ala180, Gly179, Val215, His218 | Ser172, Ala182, His183, Ala184, Glu219, His228, Arg214, Thr241, Leu235, Tyr240 |
The docking study of the mutated MMP9 showed that the 13 chosen phytochemicals possessed binding affinities between -7.7 to -9.6 kcal/mol, which showed a strong potential of interaction with the target. The control drug, doxorubicin, had a binding affinity of -8.3 kcal/mol. It was interesting to note that some of the phytochemicals exhibited a binding energy of less than and greater than that of the control, indicating a higher level of interaction and possible inhibitory activity.
Asperuloside (−7.7 kcal/mol) gave the lowest binding affinity, and Cycloartanol (−9.6 kcal/mol) gave the highest binding affinity. Table 5 is a summary of binding affinity and amino acid interactions with mutated MMP9. Figure 9b represents the 2D interaction of the compounds with the Mutant MMP9, respectively.
| Sl. No. | Compounds | Binding Affinity (Kcal/mol) | No Of H Bond | Hydrogen Bond | Carbon Hydrogen / Pi-Anion /Pi-Sulfur /Pi-Pi Stacked / Pi- Alkyl Bond | Van Der Waals |
|---|---|---|---|---|---|---|
| 1 | Beta-Amyrin | -8.4 | His411, Pro180 | Asp182, Leu187, Ala189, Ala191, Gln402, His405, His190, Gly178, Tyr179, Phe110 | ||
| 2 | Chlorogenic acid | -9.4 | 3 | Thr426, Arg424, Ala417 | His401, Val398, Leu188, Tyr423, Pro421, Pro430 | Ala189, Leu187, Gln402, Leu397, Glu427, Pro429, Gly428, Glu416, Pro415, Met419, Tyr420, Met420, His411, |
| 3 | Cycloartenol | -9.5 | 1 | Phe110 | Leu188, Val398, His401, Leu397, Tyr423 | Pro421, Gln401, Tyr420, Met422, Leu418, Arg424, Leu187, Ala189, His405, His411, His190 |
| 4 | Asperuloside | -7.7 | 3 | Tyr423, His401, His405 | Val398, Ala191 | Tyr393, Met422, Gly186, Leu188, Pro421, His411, Gln402, His190, Leu187, Ala189, Phe110, |
| 5 | Alpha-Curcumene | -7.9 | Leu418, Arg424, His401, Val398, Tyr423, Leu188 | Tyr420, Met422, Pro421, Glu416, Leu397, Gln402, Pro430, Ala417 | ||
| 6 | 15-Demethyl plumieride | -7.8 | 7 | Gly186, Met422, Tyr423, Leu188, Ala189, His405, His401 | Tyr420 | Tyr393, Leu187, Pro421, Leu418, Leu397, Arg424, Val398, Gln402, His411 |
| 7 | Delta7-Avenasterol | -8.1 | Leu187, His190, Phe110, Leu188, Val398, His401 | His411, Ala189, Gln402, Arg424, Pro421, Met422, Tyr423, Leu397, Tyr420, Tyr179 | ||
| 8 | Cycloartanol | -9.6 | 1 | Phe110 | His401, Tyr423 | Pro421, His411, Leu188, Tyr420, Val398, Met422, Leu418, Arg424, Leu397, Ala189, Gln402, Leu187, His190, His405, Glu111 |
| 9 | Lanosterol | -8.6 | Tyr423, His401, Phe110 | His411, Pro421, Tyr420, Arg424, Leu418, Leu397, Val398, Met422, Leu188, Gln402, Ala189, His190, Leu187 | ||
| 10 | (3S)-3-hydroxy-2,3-dihydro-1H-pyrrolo [2,1-b]quinazolin-9-one | -9.1 | 1 | Met422 | Leu418, Leu397, Arg424, His401, Val398 | ALa417, Tyr420, Pro421, Leu188, Tyr423 |
| 11 | Campesterol | -9.2 | 1 | Phe110 | His411, His401, Leu418 | His190, Leu187, Ala189, Gln402, Val398, Arg424, Tyr420, Tyr423, Leu397, Met422, Leu188, Pro421 |
| 12 | beta-Sitosterol | -8.6 | His401, Leu187 | Leu418, Tyr420, Tyr423, Leu397, Val398, Pro421, Leu188, His411, Tyr179, Phe110, His190, His405, Ala189, Met422, Gln402 | ||
| 13 | Stigmasterol | -9 | His401, Val398, Leu188, Leu187 | His411, Pro421, Tyr423, Tyr179, Phe110, His190, Ala189, Gln402, Met422, Tyr420, Leu418, Leu397, | ||
| C | Doxorubicin | -8.3 | Pro421 | Phe110, His411, His405, Gln402, His401, Tyr420, Met422, Val398, Tyr423, Leu188, Leu187, Gly186, Tyr393 |
In docking studies, it was established that all 13 phytochemicals showed good binding affinities with the mutated MMP1 and MMP9 proteins, meaning that they could be effective modulators of those targets.
Analysis of dynamic simulations
Molecular Dynamics (MD) simulations were conducted to evaluate the structural stability, conformational dynamics, and binding interactions of Delta-7-Avenasterol with both wild-type and mutated MMP-1 and MMP-9 proteins over a 100 ns trajectory.
Root Mean Square Fluctuation (RMSF)
The RMSF analysis demonstrated that the MMP9 wild type maintained most residue fluctuations below 0.3 nm, with moderate peaks up to 0.5 nm observed in terminal and loop regions, indicating overall structural stability with expected flexibility at disordered sites. In contrast, the MMP9 mutant displayed markedly higher structural mobility, with pronounced spikes exceeding 1.0 nm, particularly between atoms 100-200 and around 400, suggesting mutation-induced local flexibility and adaptive motions in regions proximal to the binding cleft. This enhanced mobility may promote ligand accommodation but could also reduce long-term interaction persistence (Figure 10a).
Protein and Ligand RMSD
Protein backbone RMSD analysis revealed that MMP9 wild type exhibited the lowest displacement (0.15-0.20 nm post-equilibration), reflecting a highly stable complex. The MMP9 mutated system displayed slightly higher RMSD values (0.20-0.35 nm), indicative of greater conformational movement induced by mutation. Both MMP1 variants showed intermediate stability, with RMSD values fluctuating between 0.20-0.30 nm. Ligand RMSD trends further supported these findings: Delta-7-Avenasterol remained tightly anchored within the MMP9 wild-type binding pocket (0.15-0.25 nm), confirming persistent interaction. MMP9 and MMP1 mutants showed slightly higher ligand RMSD values (0.10-0.20 nm, with occasional spikes), while the MMP1 wild type demonstrated the greatest ligand mobility, with transient excursions approaching 0.25 nm. These results suggest that wild-type proteins generally confer more stable ligand retention, whereas mutations introduce flexibility and reduced binding constriction (Figure 10b, c & d).
Radius of gyration (Rg)
The Rg profiles indicated steady structural compactness across all systems, with subtle differences between wild-type and mutated forms. MMP9 wild type exhibited the highest and most stable Rg (~1.65 nm), implying a more extended yet stable conformation. In comparison, MMP1 and mutated complexes showed slightly lower Rg values (1.50-1.55 nm), consistent with tighter overall packing (Figure 10e).
Solvent Accessible Surface Area (SASA)
Consistent with Rg behavior, the SASA analysis revealed that MMP9 wild type maintained the largest solvent-exposed surface area (110-115 nm²) throughout the simulation, suggesting a more open and flexible conformation. In contrast, the mutated and MMP1 complexes exhibited smaller SASA ranges (85-100 nm²), reflecting comparatively compact and less solvent-accessible structures (Figure 10f).
Hydrogen Bond Analysis
Hydrogen bond evaluations showed that MMP9 wild type sustained 1-2 hydrogen bonds consistently across the trajectory, ensuring stable ligand engagement. Mutated forms, particularly MMP9, formed fewer or less persistent hydrogen bonds, indicating mutation-driven alterations in interaction networks that modestly decreased overall binding stability (Figure 10 g).
Per-Residue Interaction Mapping
Key interacting residues contributing to Delta-7-Avenasterol stabilization were identified across all complexes. In the MMP1 wild type, residues A182, H218, Y240, and S172 established strong hydrogen bonding and hydrophobic contacts. The MMP1 mutated complex primarily relied on V215, Y237, H222, and P238 as stabilizing residues, reflecting adaptive reorganization of the binding environment. For MMP9 wild type, the main interacting residues included H150, D148, G117, T102, and R119, whereas MMP9 mutated engaged H401, Q402, Y423, L418, and V398, indicating a shifted binding interface driven by mutation.
Overall Stability Assessment
Collectively, the MD simulations confirmed that Delta-7-Avenasterol maintains the most stable binding interactions with MMP9 wild type, supported by low RMSD values, compact yet flexible topology (Rg), and persistent hydrogen bonding. A mutation in MMP9 introduced localized flexibility, while MMP1 variants displayed moderate dynamic stability, with the MMP1 wild type being the most fluctuating. These findings highlight the structural resilience of Delta-7-Avenasterol within MMP9 wild type, suggesting this target as the most favorable binding environment among the evaluated systems.
RT-PCR Analysis of MMP1 and MMP9 Expression
The transcriptional response of matrix metalloproteinases MMP1 and MMP9 was examined in KB oral carcinoma cells following 24 hr treatment withG. aparine at its IC₅₀ concentration (272.1 µg/mL). Treatment withG. aparine resulted in a pronounced reduction in MMP1 mRNA expression, with treated cells exhibiting a relative expression level of 0.41 compared to control cells. Similarly, MMP9 expression was decreased following treatment, with transcript levels reduced to 0.59 relative to untreated cells. The concurrent down-regulation of MMP1 and MMP9 indicates suppression of genes associated with extracellular matrix degradation and invasive potential in KB oral carcinoma cells (Figure 11).
DISCUSSION
Homeopathy has long been implemented as a supplemental treatment option for chronic conditions, including cancer (Frenkel, 2015). Though its mechanism of action is still debatable, a few studies have demonstrated that it could be beneficial in curing certain illness symptoms while improving quality of life. A systematic review by Wagenknecht et al., analysed 18 studies including more than 2000 patients with multiple forms of cancer, like gastrointestinal cancer, breast tumours, lung tumours, haematological cancer, head as well as neck tumours, renal cell cancer, pancreatic cancer, etc., (Wagenknecht et al., 2023). The patients were administered single or combination homeopathic preparations systemically or as mouth rinses. The research examined the outcomes, including therapy-induced toxicity, functional well-being, life expectancy, mental and emotional health, and tolerability. The outcomes varied, with few studies demonstrating the advantages, while others reported no significant improvement. Another study, by Frass et al., involving 410 cancer-affected patients undergoing integrative homeopathic treatment in addition to conventional therapy, demonstrated notably enhanced health condition and quality of life relative to a control group (Frass et al., 2015). This investigation suggests the prospective effects of homeopathic treatment in improving well-being alongside an established cancer regimen. In spite of insufficient therapeutic substantiation, preliminary studies suggest possible biological activity in high-dilution homeopathic formulations that could include modulation of apoptotic pathways in cancer systems (Dos Santos et al., 2021). Homeopathy treatment has been examined in the context of cancer care, where the treatment has been found to have a complementary effect in conventional cancer treatment. Certain formulations showed positive response in fatigue, anxiety, and tolerance to therapy. This fact suggests that, despite not being an alternative method to traditional therapy, the targeted homeopathic formulations may provide complementary pain management in managing pain and increasing compliance with treatment in the context of the conventional oncology therapy. The clinical findings have demonstrated the quality of life and reduced treatment side effects in people who receive supportive therapies of homeopathy compared to those who receive treatment using conventional therapy alone (Rostock et al., 2011).
G. aparine L., commonly known as Cleavers or catchweed, is an annual herb belonging to the family Rubiaceae and has traditionally been used in herbal and alternative medicine as a skin condition, lymphatic, and detoxifying agent. Homeopathically, it is regularly prescribed in the case of lymph node enlargement, urinary tract issues, and long-term dermatology issues, but empirical studies to support the assertion are scarce. The recent phytochemical studies have revealed the presence of a complex range of iridoids, flavonoids, phenolic acids, and lipophilic compounds that prove its therapeutic effects (Goryacha et al., 2014; Ilina et al., 2019; Shynkovenko et al., 2018). The pharmacological use of G. aparine has been expanded to include anti-oxidant, anti-microbial, anti-tumour, and immunostimulant activity in vitro and in vivo (Ilina et al., 2020). It was shown to have a significant wound healing, antimicrobial, and antioxidant effect, which confirmed its use in ethnomedicine (Beirami et al., 2024). G. aparine extracts have been shown to induce programmed cell death, repress the expression of angiogenic cytokines such as VEGF and IL-8, and disrupt ERK1/2 signalling of breast cancer and leukaemia cell lines, hence demonstrating the anticancer potential of this herb (Atmaca, 2017; Shi et al., 2016). Yoon et al. (2005) established that crude extracts of G. aparine inhibited metastatic development to lungs linked with stimulation of immune response, yet lacked direct cytotoxic impact on B16-BL6 melanoma cells (Yoon et al., 2005). Overall, the evidence underpins the renewed interest in G. aparine as a candidate in holistic cancer treatment and phytomedicine, which makes it necessary to conduct further experimental and patient-related studies.
In the present study, phytochemicals from G. aparine were analysed for their anticancer potential against OSCC using an in silico approach. This integrative computational approach revealed several significant patterns and biological insights. Initially, 38 phytochemicals were curated using the IMPPAT database (Mohanraj et al., 2018; Vivek-Ananth et al., 2023).
Drug-like properties of 38 phytochemicals of 344 were evaluated based on the Rule of Five. The majority of compounds met important pharmacokinetic and physicochemical requirements, which suggested favourable drug-like properties and oral bioavailability. These findings reinforce the therapeutic importance of the discovered phytochemicals (Lipinski et al., 2001; Lipinski et al., 2012). Pharmacokinetic studies as a predictive of blood-brain barrier permeability, gastrointestinal absorption, P glycoprotein substrate behaviour and cytochrome P450 isozyme inhibition, offered a clue into the safety and systemic bioavailability of compounds (Singh, 2006). End toxicity testing using ProTox III identified the possibility of a carcinogenic, mutagenic, hepatotoxic, and cytotoxic liability, and narrowed down a list of 13 phytochemicals that had no anticipated toxicity phenotype (Banerjee et al., 2018).
The selection of target genes was based on publicly available datasets of gene-disease associations (GeneCards, DisGeNet, GEO) that provided an intersection of 17 genes that are highly involved in OSCC pathogenesis. MMP1 and MMP9 were then determined as core nodes in the network pharmacology assessment using STRING and Cytoscape, using centrality measures like network degree. These findings are consistent with the large body of literature recognising MMP1 and MMP9 as key mediators of extracellular matrix degradation, tumour invasion, angiogenesis, and metastatic progression, and their dual targeting in OSCC. These enzymes are essential in the reorganisation of the extracellular matrix, tumour invasion and metastasis, and therefore they are important targets in cancer treatment. It is indicated that MMP1 and MMP9 are upregulated in the tissues and saliva samples of OSCC, and it is assumed that they serve as diagnostic markers and a therapeutic target (Cai et al., 2022; Chang et al., 2020; Hema Shree et al., 2019).
Expression profiling further corroborated the relevance of these targets: both genes exhibited marked upregulation in tumour tissues compared with normal controls, with a strongly positive co-expression trend suggesting possible co-regulation or shared regulatory pathways. This is also consistent with the GEPIA data reported from a study by Byju et al., showing high expression of MMP1, MMP2, and MMP14 in HNSC, followed by a few other MMPs, including MMP9 (Byju et al., 2024). The observed positive association between MMP1 and MMP9 expression indicates potential co-regulation and interaction in cancer progression. In terms of prognosis, elevated MMP1 and statistically insignificant MMP9 expression correlated with poor overall survival with unfavourable prognosis in HNSC. Reinforcing their role as prognostic indicators and drug targets, comparable prognostic patterns for MMP1 and MMP9 have been documented in breast cancer and uveal melanoma cohorts (Cheng et al., 2022; Wang et al., 2021). Likewise, expression of MMP1 and MMP9 was positively correlated with immune cell infiltration CD4+.
Molecular docking studies predicted that all 13 phytochemicals possess strong binding affinity toward MMP1, with 4 of them also showing significant affinity for MMP9 in OSCC. MMP1 and MMP9 are two key therapeutic targets known to play significant roles in the progression of OSCC. Genetic alterations in the target genes MMP1 and MMP9 were identified using cBioPortal. Mutations have been detected in MMP1, indicating its potential involvement in tumour progression and genetic instability in OSCC. In contrast, MMP9 exhibits comparatively lower mutation frequency but remains biologically significant due to its functional role in extracellular matrix degradation and metastasis. To evaluate the potential of the phytochemicals to target these mutant forms, molecular docking was also performed between the mutated MMP1 and MMP9 proteins and the 13 selected phytochemicals (Gao et al., 2013, 2022).
A detailed analysis was conducted using Swiss Similarity to evaluate the structural similarity between lead and known drugs (Gfeller et al., 2014). The lead candidates exhibited high structural homology to FDA-authorized drugs, highlighting the translational value of these natural compounds.
When the phytochemical binding efficiency was studied against the wild and mutated gene targets, it was noted that most of the compounds docked with the mutated MMP1 had binding affinities similar to the high-affinity compounds with the wild-type target. Conversely, a number of compounds showed significantly better binding affinities to mutated MMP9, which implies a higher potential to interact as compared to their wild-type counterparts, and indicates a higher therapeutic potential. Eleven compounds had good affinity with the wild-type and mutant MMP1, and one compound had good affinity with both the wild-type and mutant MMP9. Asperuloside was identified to have high binding affinities with the wild type and the mutant MMP1 and MMP9. Molecular dynamics simulation also points to the fact that Delta-7 Avenesterol binding is stabilized by the coordinated network of hydrogen bonds and hydrophobic contacts, which includes extremely conserved amino acids of the wild-type proteins. The wild-type forms are described by the lower values of RMSD and stable values of RMSF, which are compatible with the cohesive interaction landscapes that facilitate the high-affinity binding of the ligands and low conformational change. Conversely, mutants bring in local flexibility (as seen in marked RMSF peaks and increased SASA), weakening hydrogen bond stability, and increasing variability in the dynamics of the ligand, as evidenced by higher RMSD spikes and less compact protein folds.
The major changes in the RMSF profiles of mutant proteins and, in particular, localized flexibility in the loop region that surrounds the active site contribute to the explanation of the observed variability in the ligand accommodation. The 2D interaction maps of 0, 50, and 100ns snapshots also validate the existence and change of key binding contacts as time progresses. This change with time is indicative of the dynamic plasticity of the pocket and ligand, which is regulated by structural mutation. These results are consistent with previous computational investigations of metalloproteinases, which have shown that the flexibility of a pocket and rearrangements caused by mutations can greatly influence the ligand selectivity and affinity.
Primarily, Solvent Exposure (SASA) and protein compactness (Radius of gyration) patterns of MMP9-WT and mutant and MMP1 complexes highlight how an increased structural openness can allow easy ligand access but not necessarily stabilize the interaction. The decreased hydrogen bonding of mutants highlights the fact that non-polar and flexible interactions assume compensatory functions in protein-ligand recognition in regions that are normally recognized by classical, conserved hydrogen bonds.
Besides the computational analyses, transcriptional response of the most significant extracellular matrix-remodelling genes, MMP1 and MMP9, in KB oral carcinoma cells after the treatment with G. aparine was studied. The RT-PCR showed that the expression of MMP1 mRNA was significantly reduced after 24 hr of treatment, and MMP9 mRNA was also reduced in comparison with untreated control cells. This is consistent with previous results that extracts of several East Asian botanical extracts can transcriptionally suppress MMP1, and with studies that have shown extracts of Salvia miltiorrhiza to prevent MMP9 expression and invasion of cancer cells (Kim et al., 2017; Kim et al., 2007). Remarkably, these transcriptional alterations were achieved in the conditions of the experiment, which also resulted in substantial cytotoxic effects on KB cells, as can be seen in the MTT assay. Combined, these results suggest that G. aparine can modulate OSCC-relevant transcriptional-level molecular pathways that regulate tumour invasion and matrix degradation.
The presence of such bioactive constituents in G. aparine suggests that its traditional use in managing various ailments, including inflammatory and lymphatic disorders, may be partly attributed to its phytochemical profile. The observed anticancer potential is likely a result of synergistic interactions among these compounds targeting MMP1 and MMP9. Furthermore, the biological activity of compounds was validated utilizing the PASS web server (Druzhilovskiy et al., 2017). The PASS analysis predicted high probabilities for anticancer activity, thus reinforcing the docking findings and supporting the therapeutic promise of these compounds.
The integration of molecular docking, drug-likeness screening, ADME/toxicity screening, and structural similarity analysis offers a powerful multi-level computational platform to identify promising bioactive agents of G. aparine with minimum risks of causing adverse effects on the treatment of OSCC. Moreover, the presence of compounds that display a positive structural similarity to known drugs is an added assurance in the fact that the compounds are drug-like. The methodology in this way reduces the number of phytochemicals obtained through G. aparine as a potential inhibitor of MMP1 and MMP9 in OSCC.
Overall, the in silico data suggest that the anticancer effects of the homeopathy medication G. aparine could be explained by a set of phytochemicals that could regulate the pathways mediated by MMP. The results will serve as the basis of additional experimental evidence on the anticancer efficacy of homeopathy drug - G. aparine, especially in the treatment of OSCC. In general, these results provide a good platform upon which downstream validation, including in vitro/in vivo studies, can be done and point to a time-saving pathway of natural compound libraries to therapeutically relevant OSCC agents.
CONCLUSION
This paper provides a broad in silico analysis of the G. aparine-derived phytochemicals with respect to their therapeutic implications in attacking MMP1 and MMP9, which are two key genes in OSCC. It is worth noting that 13 bioactive compounds were discovered, which have positive pharmacological properties, low toxicity, and acceptable drug-like properties, as future safe and naturally produced therapeutic agents. Moreover, our research hypothesized that the bioactive compound Asperuloside, which is a homeopathic medicine of G. aparine, can possibly be able to interact with wild-type MMP1 and MMP9 with high binding affinity and mutated MMP1 and MMP9.
Notably, G. aparine has a historical role in homeopathy and traditional medicine systems, where it has served in its detoxifying and anti-inflammatory effects. The current results combine the classical homeopathic knowledge with the contemporary computational biology, providing the scientific evidence of the anticancer capabilities of the latter at the molecular scale. This integration of the ancient homeopathic knowledge with the current drug discovery methods is leading to new avenues of the development of evidence-based alternative medicines. These findings prompt additional preclinical studies and clinical investigation of homeopathy-inspired preparations, especially those that include well-defined phytoconstituents with established molecular interactions. Future research can also be done to test the synergy of such compounds with traditional therapies to improve the effects of treatment and decrease the toxicity of such compounds in OSCC and other malignancies.