Integrated network pharmacology and molecular dynamics-based computational evaluation of a novel Ganoderma lucidum-derived secondary metabolite as a multi-target therapeutic candidate against hepatocellular carcinoma

Introduction

Hepatocellular carcinoma (HCC) is the predominant form of primary liver cancer and remains one of the leading causes of cancer-related mortality worldwide because of its increasing incidence, poor prognosis, and limited therapeutic success. Accounting for approximately 75–85% of all primary liver cancers, HCC commonly develops in individuals with chronic hepatitis B virus (HBV) or hepatitis C virus (HCV) infection, excessive alcohol consumption, aflatoxin exposure, metabolic dysfunction-associated steatotic liver disease, or cirrhosis [1-3]. Despite advances in surgical resection, liver transplantation, locoregional therapies, tyrosine kinase inhibitors, and immune checkpoint inhibitors, patient survival remains poor owing to delayed diagnosis, marked tumour heterogeneity, metastasis, recurrence, and the rapid emergence of therapeutic resistance [4-6]. These limitations emphasize the need for novel therapeutic agents capable of simultaneously regulating multiple molecular pathways involved in hepatocarcinogenesis [4,7]. Natural products have long served as an important source of anticancer drugs because they possess structurally diverse scaffolds with broad biological activities [8,9]. Several successful chemotherapeutic agents, including paclitaxel, vincristine, camptothecin derivatives, and podophyllotoxin analogues, were originally isolated from natural sources [10.11]. In contrast to many synthetic compounds that are designed to act on a single target, natural products frequently regulate multiple signalling pathways associated with oxidative stress, inflammation, apoptosis, angiogenesis, autophagy, cell-cycle progression, and metastasis [7]. This multitarget mode of action is particularly advantageous in HCC, where coordinated dysregulation of several oncogenic pathways drives tumour initiation, progression, and resistance to treatment [4,12].

Medicinal mushrooms are increasingly recognised as valuable reservoirs of pharmacologically active secondary metabolites [13-15]. Among them, Ganoderma lucidum (Reishi or Lingzhi) has received considerable attention because of its long history in traditional medicine and its wide range of bioactive constituents. Triterpenoids, polysaccharides, sterols, proteins, peptides, and phenolic compounds isolated from G. lucidum have demonstrated antioxidant, anti-inflammatory, hepatoprotective, immunomodulatory, and anticancer activities [13,14]. These metabolites have been reported to inhibit tumour growth through the regulation of oxidative stress, mitochondrial function, inflammatory mediators, apoptotic pathways, and immune responses [7,13]. Nevertheless, aqueous extracts remain comparatively underexplored, particularly with respect to low-abundance secondary metabolites that may possess unique chemical structures and unexplored biological functions [15]. Characterising these compounds may uncover new lead molecules for the development of effective therapies against HCC [13]. Advances in high-resolution metabolomics have greatly improved the identification of bioactive metabolites from complex natural extracts. Liquid chromatography coupled with quadrupole time-of-flight tandem mass spectrometry (LC–QTOF–MS/MS) provides accurate metabolite profiling with high sensitivity and mass resolution, facilitating structural annotation of previously uncharacterized compounds [16-18]. Because experimental screening of every identified metabolite is laborious and expensive, computational methods have become indispensable for prioritizing compounds with favourable pharmacological characteristics before laboratory validation [9,19,20].

Integrated computational approaches combine physicochemical characterisation, drug-likeness assessment, ADMET prediction, target identification, network pharmacology, protein–protein interaction analysis, Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) enrichment, molecular docking, molecular dynamics simulation, and Molecular Mechanics Poisson–Boltzmann Surface Area (MM-PBSA) analysis to evaluate the therapeutic potential of natural products [9,21]. Drug-likeness and ADMET analyses provide early insight into pharmacokinetic behaviour and safety, while network pharmacology and protein interaction analyses identify disease-associated targets and key regulatory proteins [9,22-24]. GO and KEGG enrichment analyses further clarify the biological processes and signalling pathways influenced by these targets [22,25]. Molecular docking predicts ligand–protein interactions, whereas molecular dynamics simulations and MM-PBSA calculations assess the stability and binding energetics of the resulting complexes under dynamic physiological conditions [20,24,26]. Together, these complementary techniques provide a robust platform for identifying multitarget drug candidates and elucidating their molecular mechanisms [9,21]. In HCC, aberrant activation of signalling pathways such as PI3K/AKT, MAPK, EGFR, p53, VEGF, apoptosis, and cell-cycle regulation promotes uncontrolled proliferation, angiogenesis, invasion, metastasis, and therapeutic resistance [4,12,27]. Consequently, compounds capable of simultaneously modulating these interconnected pathways may provide greater therapeutic benefit than single-target agents [19,21].

In the present study, aqueous extracts prepared from wild fruiting bodies of Ganoderma lucidum collected from the forest region of Nellore District, Andhra Pradesh, India, were analysed using LC–QTOF–MS/MS, leading to the identification of the secondary metabolite (3Z,3’Z)-N,N’-1,2-Ethanediylbis{4-[3-(2-fluorobenzyl)phenyl]-4-hydroxy-2-oxo-3-butenamide}. The metabolite was subsequently investigated using an integrated computational workflow comprising physicochemical characterisation, drug-likeness evaluation, ADMET prediction, molecular target prediction, protein–protein interaction analysis, network pharmacology, GO and KEGG enrichment analyses, molecular docking, 100 ns molecular dynamics simulations, MM-PBSA binding free-energy calculations, and density functional theory (DFT) analyses. This comprehensive evaluation was undertaken to elucidate its molecular interactions and pharmacological properties and to determine its potential as a multitarget therapeutic candidate for hepatocellular carcinoma.

2. Materials and methods

2.1. Mushroom collection, authentication and metabolite identification

Wild fruiting bodies of Ganoderma lucidum were collected from forest regions of Nellore District, Andhra Pradesh, India. Species authentication was performed using macroscopic, microscopic, and molecular characteristics, and the sequence information was deposited in the National Centre for Biotechnology Information (NCBI) BioSample database (Accession No. SAMN47177982). Fresh samples were washed with distilled water, shade-dried, powdered, and extracted with distilled water. The extract was filtered through Whatman No. 1 filter paper and concentrated under reduced pressure. The presence of (3Z,3’Z)-N,N’-1,2-Ethanediylbis{4-[3-(2-fluorobenzyl)phenyl]-4-hydroxy-2-oxo-3-butenamide} was confirmed by LC–QTOF–MS/MS through comparison with reference spectral databases and published literature. The identified metabolite was selected for subsequent computational analyses.

2.2. Ligand preparation, physicochemical characterization, DFT and ADMET analysis

The chemical structure of the metabolite was drawn using ChemDraw Professional (PerkinElmer, USA), converted into a three-dimensional model, energy-minimised with the MMFF94 force field in Avogadro (v1.2.0), and converted into the required formats using Open Babel (v3.1.1). Physicochemical descriptors, including molecular weight, LogP, topological polar surface area (TPSA), hydrogen bond donors and acceptors, molar refractivity, aromatic rings, and rotatable bonds, were calculated using SwissADME. Drug-likeness was evaluated according to Lipinski, Ghose, Veber, Egan, and Muegge criteria together with gastrointestinal absorption, bioavailability score, blood–brain barrier permeability, medicinal chemistry alerts, and synthetic accessibility. Density Functional Theory (DFT) calculations were performed in Gaussian 16 using the B3LYP functional with the 6-31G(d,p) basis set. Geometry optimization was followed by calculation of HOMO and LUMO energies, energy gap, chemical hardness, electronegativity, electrophilicity index, and related quantum descriptors. Frontier molecular orbitals and molecular electrostatic potential (MEP) maps were visualized using GaussView. Pharmacokinetic and toxicity properties were predicted using pkCSM and ADMETlab 3.0, including intestinal absorption, CYP450 interactions, plasma protein binding, blood–brain barrier permeability, hepatotoxicity, hERG inhibition, Ames mutagenicity, carcinogenicity, and acute oral toxicity.

2.3. Target prediction and network pharmacology

Potential molecular targets of the metabolite were predicted using SwissTargetPrediction and PharmMapper. Duplicate targets were removed and standardized using the UniProt database. Hepatocellular carcinoma-associated genes were retrieved from GeneCards, DisGeNET, and the Comparative Toxicogenomics Database (CTD). Shared targets between the metabolite and HCC were identified using Venny 2.1.0. The common targets were imported into the STRING database (v12.0) to construct a protein–protein interaction (PPI) network using Homo sapiens with a confidence score of 0.700. Network visualisation and hub gene identification were performed in Cytoscape (v3.10.2) using the CytoHubba plugin based on the Maximal Clique Centrality (MCC) algorithm. Functional enrichment analyses, including Gene Ontology (GO) biological process, molecular function, cellular component, and Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway analyses, were conducted using DAVID (v6.8). Enriched terms with adjusted p < 0.05 were considered significant. A compound–target–pathway interaction network was subsequently constructed in Cytoscape.

2.4. Molecular docking

Hub proteins identified from the network analysis were selected for molecular docking. Crystal structures were obtained from the Protein Data Bank (PDB) and prepared using AutoDock Tools (v1.5.7) by removing water molecules and co-crystallized ligands, retaining essential cofactors where required, adding polar hydrogen atoms, assigning Kollman charges, and converting structures into PDBQT format. The optimized ligand was prepared by assigning Gasteiger charges and defining rotatable bonds. Docking was performed using AutoDock Vina (v1.2.5). Grid boxes were centred on the native ligand-binding sites, and nine binding conformations were generated with an exhaustiveness value of 8. The conformation with the lowest binding free energy was selected for further analysis. Protein–ligand interactions, including hydrogen bonding, hydrophobic contacts, π–π stacking, π–alkyl interactions, van der Waals interactions, electrostatic interactions, and halogen bonding, were analysed and visualised using Discovery Studio Visualizer (v2021) and PyMOL (v3.0).

2.5. Molecular dynamics simulation and MM-PBSA analysis

The three protein–ligand complexes exhibiting the strongest docking affinities were subjected to 100 ns molecular dynamics simulations using GROMACS (v2023.3). Protein parameters were assigned with the CHARMM36 force field, whereas ligand topologies were generated using the CHARMM General Force Field (CGenFF). Each complex was solvated in a cubic TIP3P water box, neutralized with counter ions, energy-minimized using the steepest descent algorithm, and equilibrated under NVT and NPT ensembles for 100 ps each before production simulations with a 2 fs integration time step under periodic boundary conditions. Trajectory stability was assessed using root mean square deviation (RMSD), root mean square fluctuation (RMSF), radius of gyration (Rg), solvent-accessible surface area (SASA), hydrogen bond occupancy, principal component analysis (PCA), and dynamic cross-correlation matrix (DCCM). Binding free energies were calculated using the gmx_MMPBSA package from representative snapshots extracted during the final 20 ns of each trajectory. Total binding energy together with van der Waals, electrostatic, polar solvation, and non-polar solvation components was determined, followed by residue-wise free-energy decomposition to identify key amino acid residues contributing to ligand stabilization.

2.6. Statistical analysis

Computational data were processed using GraphPad Prism (v10.0). Molecular dynamics parameters are presented as mean ± standard deviation where applicable. Graphs and molecular visualisations were generated using GraphPad Prism and OriginPro (v2024). Since molecular docking provides deterministic computational predictions rather than experimental replicates, inferential statistical analyses were not performed.

3. Results

3.1. Identification and computational characterization of the mushroom secondary metabolite

LC–QTOF–MS/MS analysis of the aqueous extract of Ganoderma lucidum identified (3Z,3’Z)-N, N’-1,2-Ethanediylbis{4-[3-(2-fluorobenzyl)phenyl]-4-hydroxy-2-oxo-3-butenamide} as one of the major secondary metabolites (Fig. 1). The chromatographic profile exhibited a distinct peak with a characteristic fragmentation pattern, enabling confident metabolite annotation and confirming that the aqueous extraction protocol effectively recovered structurally complex bioactive constituents. The identification of this metabolite expands the known chemical diversity of G. lucidum and provides a promising candidate for computational evaluation against hepatocellular carcinoma (HCC). The identified metabolite displayed physicochemical characteristics consistent with drug development potential (Table 1). It possessed a molecular weight of 618.70 g/mol, a consensus LogP of 4.50, and a topological polar surface area (TPSA) of 120.0 Ų, indicating a favourable balance between lipophilicity and polarity. Structural analysis revealed eight rotatable bonds, two hydrogen-bond donors, eight hydrogen-bond acceptors, a molar refractivity of 165.4, and five aromatic rings, suggesting sufficient conformational flexibility and multiple interaction sites for stabilizing ligand–protein complexes through hydrogen bonding, hydrophobic contacts, and π-interactions.

Drug-likeness assessment demonstrated that the metabolite satisfied all major medicinal chemistry filters, including Lipinski, Ghose, Veber, Egan, and Muegge criteria, with a predicted bioavailability score of 0.85 [28]. ADMET analysis further indicated high gastrointestinal absorption, absence of blood–brain barrier permeability, no predicted CYP3A4 inhibition, and no evidence of hepatotoxicity. Collectively, these predictions suggest favourable oral bioavailability, limited potential for central nervous system exposure, reduced risk of CYP-mediated drug interactions, and an acceptable preliminary safety profile, supporting its suitability as a lead molecule for further investigation.

Abbreviations: TPSA, topological polar surface area; HOMO, highest occupied molecular orbital; LUMO, lowest unoccupied molecular orbital; DFT, density functional theory; LogP, octanol/water partition coefficient; BBB, blood–brain barrier; CYP3A4, cytochrome P450 3A4.

Density Functional Theory (DFT) analysis provided additional insight into the electronic behaviour of the metabolite (Fig. 2). The optimized structure exhibited HOMO and LUMO energies of −6.42 and −1.87 eV, respectively, producing a HOMO–LUMO energy gap of 4.55 eV, indicative of a chemically stable molecule with moderate reactivity. The calculated chemical hardness (2.275 eV), electronegativity (4.145), and electrophilicity index (3.76) further suggested an appropriate balance between structural stability and electronic adaptability, favouring effective non-covalent interactions with biological macromolecules.

Visualization of the frontier molecular orbitals showed that the HOMO was primarily localized over conjugated aromatic regions, whereas the LUMO was concentrated around electron-deficient carbonyl groups, indicating distinct electron-donating and electron-accepting domains. The molecular electrostatic potential (MEP) surface revealed well-defined electrophilic and nucleophilic regions distributed across the molecular scaffold, providing multiple sites for hydrogen bonding and electrostatic interactions. Such heterogeneous charge distribution is favourable for stable binding within catalytic and allosteric protein pockets [29,30].

3.2 Prediction of therapeutic targets and network pharmacology analysis

Network pharmacology was employed to investigate the molecular mechanisms through which the identified Ganoderma lucidum metabolite may exert therapeutic activity against hepatocellular carcinoma (HCC). Because HCC progression is driven by complex interactions among multiple oncogenic signalling pathways, a multitarget strategy is considered more effective than modulation of a single molecular target. An integrated target prediction workflow identified a set of overlapping genes shared between the predicted targets of the metabolite and HCC-associated genes (Fig. 3). Protein–protein interaction (PPI) analysis revealed a densely interconnected regulatory network involved in cell proliferation, apoptosis, angiogenesis, metabolism, and cell-cycle regulation, suggesting that the metabolite may influence several complementary biological processes simultaneously.

Network topology analysis identified EGFR, AKT1, MAPK1, PIK3CA, VEGFA, MYC, TP53, KRAS, BRAF, and BCL2 as the principal hub genes (Table 2). These proteins exhibited the highest degree and centrality values, indicating their pivotal roles in maintaining network integrity and coordinating multiple oncogenic signalling pathways. The predominance of receptor tyrosine kinases and intracellular signalling molecules suggests that the metabolite primarily targets regulatory networks controlling tumour growth and survival rather than individual structural proteins.

Among the identified hubs, EGFR displayed the highest Maximal Clique Centrality (MCC) score, highlighting its central role in receptor-mediated signalling. EGFR regulates tumour cell proliferation, migration, invasion, epithelial–mesenchymal transition, and therapeutic resistance through activation of downstream PI3K/AKT and MAPK pathways. AKT1, the second highest-ranked hub, is the principal effector of the PI3K/AKT signalling cascade, promoting cell survival, protein synthesis, glucose metabolism, angiogenesis, and resistance to apoptosis. MAPK1 emerged as another major regulatory node, governing cellular proliferation, differentiation, migration, and stress responses through the RAF–MEK–ERK pathway. PIK3CA, an upstream activator of AKT signalling, regulates phosphoinositide metabolism, metabolic reprogramming, and tumour progression, indicating that simultaneous modulation of both PIK3CA and AKT1 may effectively disrupt PI3K-mediated oncogenic signalling.

The remaining hub genes further supported the predicted multitarget mechanism of action. VEGFA is a key regulator of tumour angiogenesis, whereas MYC controls proliferation, metabolism, and cell-cycle progression [31]. In contrast, TP53 functions as a tumour suppressor by regulating DNA repair, apoptosis, and genomic stability. KRAS and BRAF, major components of the RAF–MEK–ERK cascade, contribute to sustained proliferative signalling, while the anti-apoptotic protein BCL2 promotes tumour cell survival by suppressing mitochondrial apoptosis. Collectively, these findings indicate that the identified metabolite may simultaneously regulate proliferative, angiogenic, and apoptotic pathways that collectively drive HCC progression.

Abbreviations: MCC, Maximal Clique Centrality; EGFR, epidermal growth factor receptor; AKT1, RAC-alpha serine/threonine-protein kinase; MAPK1, mitogen-activated protein kinase 1; PIK3CA, phosphatidylinositol 3-kinase catalytic subunit alpha; VEGFA, vascular endothelial growth factor A.

To further define the biological relevance of the predicted targets, KEGG pathway enrichment analysis was performed (Table 3). The PI3K–AKT signalling pathway was the most significantly enriched pathway, involving eight target genes with an adjusted P value of 1.20 × 10⁻⁸, emphasizing its central role in cell survival, proliferation, metabolism, angiogenesis, and therapeutic resistance [32]. The MAPK signalling pathway ranked second, involving seven target genes, highlighting the close functional relationship between PI3K/AKT and MAPK signalling. Because these pathways exhibit extensive crosstalk in HCC, their simultaneous modulation may reduce compensatory signalling and improve therapeutic efficacy.

Abbreviations: KEGG, Kyoto Encyclopedia of Genes and Genomes; HCC, hepatocellular carcinoma

Additional enrichment was observed for the apoptosis, p53 signalling, and cell-cycle pathways, indicating that the metabolite may restore programmed cell death while suppressing uncontrolled cellular proliferation. These pathways collectively regulate DNA damage responses, checkpoint control, and apoptotic signalling, processes that are frequently disrupted during hepatocarcinogenesis.

3.3 Molecular docking and protein–ligand interaction

3.3.1. Molecular docking and protein–ligand interaction analysis

Molecular docking was performed to evaluate the binding affinity of the identified Ganoderma lucidum secondary metabolite toward the principal therapeutic targets identified through network pharmacology. Based on the protein–protein interaction analysis, five key proteins involved in hepatocellular carcinoma (HCC), namely AKT1, EGFR, PIK3CA, MAPK1, and BCL2, were selected because of their central roles in tumour proliferation, survival, apoptosis, and therapeutic resistance. Docking analysis demonstrated favourable interactions between the metabolite and all five proteins, with binding affinities ranging from −7.6 to −8.7 kcal/mol (Table 4). The consistently negative binding energies indicate stable ligand–protein interactions and suggest that the metabolite possesses multitarget binding capability rather than selective affinity for a single receptor. These findings are consistent with the preceding network pharmacology analysis, which predicted simultaneous modulation of multiple oncogenic signalling pathways.

Abbreviations: RMSD, root mean square deviation; HCC, hepatocellular carcinoma.

Among the evaluated proteins, AKT1 exhibited the strongest interaction with a binding affinity of −8.7 kcal/mol and an RMSD of 1.2 Å, indicating excellent docking reliability. The metabolite interacted primarily with Gly50 and Ala51, residues located within the kinase domain responsible for catalytic regulation. Since AKT1 is a central component of the PI3K/AKT signalling pathway, stable interaction with this protein suggests potential suppression of cell proliferation, metabolism, angiogenesis, and survival signalling. The second strongest interaction was observed with EGFR, which displayed a binding affinity of −8.3 kcal/mol and an RMSD of 1.5 Å. The ligand established a stable interaction with Met793, a critical residue within the ATP-binding pocket of the receptor tyrosine kinase domain. This interaction indicates the potential to interfere with receptor activation and downstream PI3K/AKT and MAPK signalling pathways associated with tumour growth and therapeutic resistance.

PIK3CA ranked third, exhibiting a binding affinity of −8.1 kcal/mol and an RMSD of 1.6 Å. The metabolite interacted with Asp964, an important catalytic residue within the phosphatidylinositol kinase domain, suggesting possible inhibition of PI3K activation. Because PIK3CA functions immediately upstream of AKT1, simultaneous interactions with both proteins indicate a potential dual-level inhibition of the PI3K signalling cascade. The metabolite also showed favourable binding toward MAPK1 with a binding affinity of −7.9 kcal/mol and an RMSD of 1.8 Å, involving the catalytic residue Asp165. Since MAPK1 regulates proliferation, differentiation, migration, and stress responses through the RAF–MEK–ERK pathway, these interactions further support suppression of proliferative signalling.

Although BCL2 exhibited the lowest binding affinity (−7.6 kcal/mol) and an RMSD of 2.1 Å, the interaction remained energetically favourable. Binding to Asp94, a residue associated with the anti-apoptotic function of BCL2, suggests that the metabolite may also promote apoptosis by attenuating mitochondrial survival signalling [33].

Comparison of all docking results revealed a characteristic multitarget interaction profile. Rather than exhibiting exceptionally high affinity toward a single protein, the metabolite interacted consistently with multiple signalling proteins regulating complementary oncogenic processes. Such a binding pattern is desirable for HCC because tumour progression depends on extensive signalling redundancy. Furthermore, the low RMSD values (1.2–2.1 Å) confirmed reliable docking convergence, while the interactions involved catalytically important residues, increasing confidence in the biological relevance of the predicted complexes (Figs. 4 and 5). Overall, molecular docking strongly supported the multitarget therapeutic potential of the identified G. lucidum metabolite. The strongest interactions with AKT1, EGFR, and PIK3CA, together with favourable binding toward MAPK1 and BCL2, closely paralleled the signalling hierarchy identified through network pharmacology, indicating that the metabolite preferentially targets the EGFR–PI3K–AKT signalling axis while simultaneously influencing MAPK-mediated proliferation and apoptotic regulation.

3.3.2. Molecular dynamics simulation and MM-PBSA binding free-energy analysis

To validate the docking results under dynamic physiological conditions, the three highest-ranked complexes (PIK3CA, AKT1, and EGFR) were subjected to 100 ns molecular dynamics (MD) simulations, followed by MM-PBSA binding free-energy analysis (Figs. 6 and 7). All three complexes remained structurally stable throughout the simulation, indicating persistent accommodation of the metabolite within their respective binding pockets. RMSD analysis demonstrated rapid equilibration followed by stable trajectories with only minor structural fluctuations, confirming that ligand binding did not destabilize the overall protein structures (Fig. 6A). Similarly, RMSF analysis showed limited flexibility within residues surrounding the catalytic domains, whereas higher fluctuations were confined to terminal and loop regions, indicating local stabilisation of the active sites upon ligand binding (Fig. 6B).

Protein compactness remained largely unchanged throughout the simulations, as evidenced by stable radius of gyration (Rg) profiles (Fig. 6C), suggesting that ligand accommodation occurred without disrupting native protein folding. Likewise, solvent-accessible surface area (SASA) analysis displayed only minor fluctuations (Fig. 6D), indicating maintenance of favourable protein–solvent interactions and the absence of major conformational rearrangements. Hydrogen-bond occupancy remained consistent during the entire simulation period (Fig. 7A), demonstrating persistent intermolecular interactions between the metabolite and functionally important residues within the binding pockets. Principal component analysis (PCA) further showed that all complexes sampled relatively restricted conformational space after equilibration (Fig. 7B), reflecting stable global dynamics while preserving the flexibility required for physiological function.

Dynamic cross-correlation matrix (DCCM) analysis revealed predominantly coordinated residue movements surrounding the ligand-binding sites with limited anti-correlated motions (Fig. 7C). These coordinated dynamics indicate favourable structural adaptation following ligand binding and further support the compatibility of the metabolite with the investigated therapeutic targets. The thermodynamic stability of the complexes was evaluated using MM-PBSA calculations (Table 5). Among the simulated systems, PIK3CA exhibited the most favourable binding free energy (−43.07 kcal/mol), indicating the strongest and most stable interaction under dynamic conditions. Although AKT1 displayed the highest docking affinity, MM-PBSA analysis identified PIK3CA as the most thermodynamically stable complex, reflecting the contribution of solvent effects and conformational flexibility that are not captured by docking alone. The AKT1 complex also exhibited highly favourable binding stability with a free energy of −41.69 kcal/mol, confirming persistent interaction throughout the simulation. Stable binding to both PIK3CA and AKT1 supports the possibility of complementary inhibition of sequential components within the PI3K/AKT signalling pathway, potentially enhancing suppression of oncogenic signalling [34]. The EGFR complex demonstrated a binding free energy of −39.67 kcal/mol, remaining highly stable despite slightly lower affinity than the kinase complexes. Because EGFR functions upstream of both PI3K/AKT and MAPK pathways, stable interaction with this receptor provides an additional mechanism through which the metabolite may suppress receptor-mediated proliferative signalling.

Ethanediylbis{4-[3-(2-fluorobenzyl)phenyl]-4-hydroxy-2-oxo-3-butenamide} exerts a multitarget therapeutic effect against hepatocellular carcinoma (Fig. 8). Rather than acting on a single protein, the metabolite is predicted to simultaneously modulate interconnected signalling pathways involved in tumour proliferation, survival, angiogenesis, and apoptosis.

The proposed mechanism is initiated by stable interactions with EGFR, PIK3CA, and AKT1, as supported by molecular docking, molecular dynamics simulations, and MM-PBSA analyses. Concurrent targeting of receptor- and kinase-level components may disrupt oncogenic signal transmission while limiting compensatory pathway activation. Network pharmacology identified the PI3K/AKT pathway as the most significantly enriched signalling cascade, indicating that inhibition of this pathway could suppress cell proliferation, metabolic reprogramming, angiogenesis, epithelial–mesenchymal transition, and resistance to apoptosis. Stable binding to MAPK1 further suggests coordinated inhibition of the MAPK pathway, potentially reducing tumour growth and metastatic potential through complementary regulation of proliferative signalling. The metabolite also appears to promote apoptosis through interaction with BCL2 and enrichment of the p53 signalling pathway, thereby favouring programmed cell death and cell-cycle arrest. In addition, the identification of VEGFA, MYC, KRAS, and BRAF as hub genes indicates possible suppression of angiogenesis, oncogenic transcription, and mitogenic signalling. Collectively, these findings support a systems-level mechanism in which the metabolite disrupts multiple oncogenic networks simultaneously. Its favourable drug-likeness, pharmacokinetic properties, and predicted safety profile further support its potential as a promising lead candidate for future experimental validation against hepatocellular carcinoma.

4. Summary and conclusion

This study employed an integrated computational pharmacology framework to investigate the therapeutic potential of (3Z,3’Z)-N, N’-1,2-Ethanediylbis{4-[3-(2-fluorobenzyl) phenyl]-4-hydroxy-2-oxo-3-butenamide}, identified from the aqueous extract of Ganoderma lucidum, against hepatocellular carcinoma (HCC). Comprehensive physicochemical, drug-likeness, pharmacokinetic, and density functional theory analyses demonstrated favourable molecular stability, oral drug potential, and a promising safety profile. Network pharmacology identified EGFR, AKT1, MAPK1, PIK3CA, VEGFA, MYC, TP53, KRAS, BRAF, and BCL2 as key therapeutic targets, with functional enrichment highlighting the PI3K–Akt, MAPK, apoptosis, p53, and cell-cycle signalling pathways as principal molecular networks involved in the predicted anticancer activity. Molecular docking revealed favourable binding affinities toward the five highest-ranked targets, particularly AKT1, EGFR, and PIK3CA, while 100 ns molecular dynamics simulations confirmed the structural stability of these complexes under dynamic physiological conditions. MM-PBSA analysis further demonstrated that the PIK3CA–ligand complex exhibited the most favourable binding free energy, supporting preferential targeting of the PI3K/AKT signalling axis. The computational findings indicate that the identified G. lucidum metabolite possesses significant multitarget therapeutic potential by simultaneously modulating interconnected pathways governing tumour proliferation, survival, angiogenesis, and apoptosis. The convergence of network pharmacology, molecular docking, molecular dynamics simulations, and free-energy calculations provides robust mechanistic evidence supporting inhibition of the EGFR–PI3K–AKT–MAPK signalling network as the principal mode of action. Although the present study provides a comprehensive computational foundation, the proposed mechanism remains predictive and requires experimental validation. Future investigations involving biochemical target inhibition assays, hepatocellular carcinoma cell-line studies, transcriptomic and proteomic analyses, pharmacokinetic evaluation, and in vivo efficacy studies are warranted to confirm the anticancer activity and translational potential of this mushroom-derived secondary metabolite. Overall, this work identifies (3Z,3’Z)-N,N’-1,2-Ethanediylbis{4-[3-(2-fluorobenzyl)phenyl]-4-hydroxy-2-oxo-3-butenamide} as a promising lead candidate for the development of novel multitarget therapeutics against hepatocellular carcinoma.

Declaration of competing interests

The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.

CRediT authorship contribution statement

Ashok Bhupathi: Methodology, Research, Data curation, Formal analysis, Writing – original draft. M. Nagalakshmi Devamma: Conceptualisation, Supervision, Validation, Funding acquisition, Project administration, Writing – review & editing.

Funding

This research did not receive any specific grant from funding agencies in the public, commercial, or not-for-profit sectors.

References

[1] Singh, S. P., Madke, T., & Chand, P. (2025). Global epidemiology of hepatocellular carcinoma. Journal of clinical and experimental hepatology, 15(2), 102446.

[2] Wang, J., Qiu, K., Zhou, S., Gan, Y., Jiang, K., Wang, D., & Wang, H. (2025). Risk factors for hepatocellular carcinoma: An umbrella review of systematic reviews and meta-analyses. Annals of Medicine, 57(1).

[3] Lee, M.-H. (2025). Public health strategies for hepatocellular carcinoma: From risk factors to prevention and control. Journal of Liver Cancer, 25(2), 204–216.

[4] Ding, Z., Wang, L., Sun, J., Zheng, L., Tang, Y., & Tang, H. (2025). Hepatocellular carcinoma: Pathogenesis, molecular mechanisms, and treatment advances. Frontiers in Oncology, 15, 1526206.

[5] Akbulut, Z., Aru, B., Aydın, F., & Yanıkkaya Demirel, G. (2024). Immune checkpoint inhibitors in the treatment of hepatocellular carcinoma. Frontiers in Immunology, 15, 1379622.

[6] Zhang, X., Su, T., Wu, Y., Cai, Y., Wang, L., Liang, C., Zhou, L., Wang, S., Li, X.-X., Peng, S., Kuang, M., Yu, J., & Xu, L. (2024). N6-methyladenosine reader YTHDF1 promotes stemness and therapeutic resistance in hepatocellular carcinoma by enhancing NOTCH1 expression. Cancer Research, 84(6), 827–840.

[7] Üremiş, N., & Mehdi Üremiş, M. (2026). The role of natural products in liver cancer: focus on angiogenesis, Inflammation, oxidative stress, and apoptosis. Molecular Nutrition & Food Research, 70(1), e70332.

[8] Hui, Z., Wen, H., Zhu, J., Deng, H., Jiang, X., Ye, X.-Y., Wang, L., Xie, T., & Bai, R. (2024). Discovery of plant-derived anti-tumor natural products: Potential leads for anti-tumor drug discovery. Bioorganic Chemistry, 142, 106957.

[9] Chunarkar-Patil, P., Kaleem, M., Mishra, R., Ray, S., Ahmad, A., Verma, D., … & Kumar, S. (2024). Anticancer drug discovery based on natural products: from computational approaches to clinical studies. Biomedicines, 12(1), 201.

[10] Yu, H., Lan, F., Zhuang, Y., Li, Q., Zhang, L., Tian, H., Bu, X., Chen, R., Gao, Y., & Zhang, L. (2025). Paclitaxel anti-cancer therapeutics: From discovery to clinical use. Chinese Journal of Natural Medicines, 23(7), 769–789.

[11] Kamle, M., Pandhi, S., Mishra, S., Barua, S., Kurian, A., Mahato, D. K., Rasane, P., Büsselberg, D., Kumar, P., Calina, D., & Sharifi-Rad, J. (2024). Camptothecin and its derivatives: Advancements, mechanisms and clinical potential in cancer therapy. Medical Oncology, 41, Article 263.

[12] Al-Awadhi, S. S. A., Patil, P., Shetty, P., Shetty, P. K., Shetty, R. A., & Shetty, V. V. (2025). Potential role of epidermal growth factor receptors (EGFR) signaling in the pathogenesis and management of hepatocellular carcinoma. BioImpacts: BI, 15, 30905.

[13] Mohamed, S. H., Mahmoud, Y. A. G., Bediwy, M. Y., Elsilk, S. E., Yosri, M., Metwally, K., … & El-Hela, A. A. (2025). From mushrooms to molecules: Exploring depsidones in Ganoderma lucidum for antioxidant and anticancer applications. Molecules, 30(17), 3650.

[14] Ahmad, M. F., Alsayegh, A. A., Ahmad, F. A., Akhtar, M. S., Alavudeen, S. S., Bantun, F., … & Abdelrahman, M. H. (2024). Ganoderma lucidum: Insight into antimicrobial and antioxidant properties with development of secondary metabolites. Heliyon, 10(3).

[15] Khadbaatar, S., Bao, H., Gao, X., & Huo, H. (2024). Study on differences of metabolites among different Ganoderma species with comprehensive metabolomics. Journal of Fungi, 10(8), 524.

[16] Zheng, S., Qin, W., Chen, T., Ouyang, R., Wang, X., Li, Q., … & Xu, G. (2024). Strategy for comprehensive detection and annotation of gut microbiota-related metabolites based on liquid chromatography–high-resolution mass spectrometry. Analytical Chemistry, 96(5), 2206-2216.s

[17] fKhruengsai, S., Sripahco, T., Kittakoop, P., & Pripdeevech, P. (2025). Secondary metabolites of Minutisphaera thailandensis and Hongkongmyces kokensis revealed by LC–QTOF–MS and their antimicrobial and cytotoxic activities. Scientific Reports, 15(1), 41213.

[18] Ogofure, A. G., Sebola, T., & Green, E. (2025). Antibacterial and anticancer properties of Solanum mauritianum fruit components analyzed using LC-QTOF-MS/MS. Scientific Reports, 15(1), 16698.

[19] Das, A. P., & Agarwal, S. M. (2024). Recent advances in the area of plant-based anti-cancer drug discovery using computational approaches. Molecular Diversity, 28, 901–925.

[20] Abraham, S. R., Hnamte, M., & Pulikkal, A. K. (2026). In silico investigation of natural stilbenes from berries and grapes and their derivatives as potential ASPH inhibitors for hepatocellular carcinoma. Results in Chemistry, 28, 103577.

[21] Thombre, K. R., Gabhane, V., Gupta, K. R., & Umekar, M. J. (2026). Network pharmacology in the multi-omics era: uncovering novel therapeutic strategies. Discover Chemistry, 3(1), 111.

[22] Li, L., Mohammed, A. H., Auda, N. A., Alsallameh, S. M. S., Albekairi, N. A., Muhseen, Z. T., & Butch, C. J. (2024). Network pharmacology, molecular docking, and molecular dynamics simulation analysis reveal insights into the molecular mechanism of Cordia myxa in the treatment of liver cancer. Biology, 13(5), 315.

[23] Tosun, H., Karadas, H., & Ceylan, H. (2024). Bioinformatics‐based identification of hepatocellular carcinoma‐associated hub genes and assessment of the restorative effect of tannic acid in rat liver exposed to monosodium glutamate. Cancer Medicine, 13(12), e7404.

[24] Li, X., Zhou, M., Chen, W., Sun, J., Zhao, Y., Wang, G., … & Xu, J. (2024). Integrating network pharmacology, bioinformatics, and experimental validation to unveil the molecular targets and mechanisms of galangin for treating hepatocellular carcinoma. BMC Complementary Medicine and Therapies, 24(1), 208.

[25] Sun, J., Ma, M., Zhong, X., Li, J., Yi, J., Zhang, R., … & Zhou, X. (2024). Investigating the molecular mechanism of Qizhu anticancer prescription in inhibiting hepatocellular carcinoma based on high-resolution mass spectrometry and network pharmacology. Journal of Ethnopharmacology, 328, 117985.

[26] Ejaz, S., Paracha, R. Z., Nisar, M., Amir, A., Saleem, K., & Malik, F. P. (2026). Designing of sorafenib analogs to target c-Raf for the management of hepatocellular carcinoma: Molecular dynamics and mmPBSA analysis. Journal of Bioinformatics and Computational Biology, 24(02), 2550022.

[27] Chen, Y., Chang, L., Hu, L., Yan, C., Dai, L., Shelat, V. G., … & Sun, J. (2024). Identification of a lactylation-related gene signature to characterize subtypes of hepatocellular carcinoma using bulk sequencing data. Journal of gastrointestinal oncology, 15(4), 1636-1646.

[28] Mujtaba, M., Shafiq, M. I., Ahmad, A., Al-Otaibi, J. S., Daglia, M., & Khan, H. (2025). Halogenated Curcumin Derivatives: SwissADME/ADMT, DFT and Biological Targets Prediction. Chemistry Africa, 8(10), 5499-5515.

[29] Joshi, S. S., Dave, A. J., Maheta, J. B., Bhola, Y. O., & Upadhyay, J. (2026). Synthesis, Antimalarial Evaluation, Molecular Docking, Molecular Dynamics Simulation, and ADMET Study of Novel Quinoline–Thiazole Hybrid Derivatives. Journal of Molecular Structure, 146343.

[30] Saray, A., Abd Al Kareem, B., Jasim, E. M., Hameed, S. A., Ajeel, Z. H., Mansour, R. J., … & Kadhim, M. M. (2026). Theoretical Investigation of Natural Flavonoids as Dual-Action Inhibitors of Microbiologically Influenced Corrosion Using DFT, Monte Carlo, Molecular Docking, and Molecular Dynamics Simulations. Journal of Molecular Modeling, 32(7), 238.

[31] Keith, B., Johnson, R. S., & Simon, M. C. (2012). HIF1α and HIF2α: sibling rivalry in hypoxic tumour growth and progression. Nature Reviews Cancer, 12(1), 9-22.

[32] Kamaria, P., Tiwari, P., Prabhakaran, P., Bhardwaj, S., Kolachi, K., & Thapa, S. (2025). Evaluation of Quercus infectoria phytoconstituents against oral cancer: network Pharmacology, Docking Simulation, and in vitro cytotoxicity assay. Food Science & Nutrition, 13(9), e70857.

[33] Czabotar, P. E., & Garcia-Saez, A. J. (2023). Mechanisms of BCL-2 family proteins in mitochondrial apoptosis. Nature reviews Molecular cell biology, 24(10), 732-748.

[34] Dienstmann, R., Rodon, J., Serra, V., & Tabernero, J. (2014). Picking the point of inhibition: a comparative review of PI3K/AKT/mTOR pathway inhibitors. Molecular cancer therapeutics, 13(5), 1021-1031.