Staphylococcus aureus is a nonmotile, widely distributed, gram-positive coccus that serves as a significant human pathogen responsible for various infections, including skin and soft tissue infections, bacteremia, pneumonia, and several diseases caused by toxins [1]. Among its numerous extracellular proteins, S. aureus strains produce various potent toxins, including alpha hemolysin, enterotoxins, leukocidins, and exfoliative toxins, all directly linked to specific disease manifestations [2, 3]. Enterotoxins produced by S. aureus are of particular concern in food toxicology because they are highly stable and remain active in foods even after processing or cooking [4]. These enterotoxins could withstand pasteurization-like heat treatments, survive acidic conditions of fermented foods, and resist proteolytic degradation in the digestive tract [5].

Staphylococcal enterotoxin B (SEB) serves as a standard example of a potent superantigen not associated with the enterotoxin gene cluster (egc) [6, 7]. It is classified as a category B select agent due to its high potency [8]. Importantly, its role as a foodborne toxin has made it a persistent challenge for toxicologists and food safety researchers [9]. Unlike synthetic chemicals, much smaller amounts of staphylococcal enterotoxin are required to elicit a toxic response [10]. Additionally, SEB is highly stable and can be easily produced in significant quantities [11]. Even at low doses, SEB has the potential to induce multi-organ system failure and can be fatal [12].

The main targets of SEB are the T cell receptors found on T cells and the MHC class II molecules present on antigen-presenting cells (APCs) [13]. This interaction results in the formation of a ternary complex involving MHC class II molecules and specific Vβ chains of the TcR due to cross-linking [14, 15]. SEB attaches to the MHC molecule outside of the peptide-binding groove without needing prior processing, activating one of the seven Vh subclasses of the TcR [16]. The activated T cells produce significant quantities of cytokines, such as interleukin 2 (IL-2), tumor necrosis factor alpha (TNF-α), and gamma interferon (IFN-γ), which lead to hyperproliferation and eventual depletion [17-19]. The released cytokines exert strong effects, causing fever, low blood pressure, dysfunction in multiple organs, and ultimately, fatal shock [20].

Given the toxicological significance of SEB, there is an urgent need for natural food-derived inhibitors that can reduce its biological activity. Natural products have historically served as the primary source of many active pharmaceutical ingredients [21]. Flavonoids, a significant class of natural compounds, are associated with a wide range of health benefits and play a crucial role in various applications within the nutraceutical, pharmaceutical, medicinal, and cosmetic industries [22, 23]. Aromadendrin is a bioactive flavonoid extracted from Pinus sibirica, Afzelia bella, and Chioanathus retusus [24-26]. Aromadendrin exhibits radical scavenging properties, as well as anti-tumor and anti-inflammatory effects [26-28]. Additionally, aromadendrin had anti-diabetic [29], immunomodulatory influence [30], anti-oxidative properties [30], and neuroprotective effects [31].

Bergenia ligulata is a perennial herb indigenous to the Himalayan region, classified within the Saxifragaceae family. It is predominantly found in rocky terrains of North India and commonly referred to as Pashanbheda within the Indian context [32, 33]. The rhizomes of B. ligulata possess a range of pharmacological properties, including diuretic, tonic, astringent, and laxative effects [34]. The chemical constituents of this species comprise bergenin, β-sitosterol, β-sitosterol-D-glucoside, and several phenolic compounds, such as afzelechin (AZC) and catechin, along with various other compounds [35]. AZC is recognized for its antioxidant, α-glucosidase inhibitory activity, anticancer, antimicrobial, antifungal activities, and its cardiovascular protective properties [36-40].

By employing molecular docking and molecular dynamics simulations, this study bridges food science and toxicology to evaluate the structural interactions between SEB and these two flavonoids. Understanding these interactions not only provides insights into toxin inhibition at the molecular level but also supports the development of food-safe, natural inhibitors that can be incorporated into dietary strategies or food preservation systems to enhance public health safety.

MD was performed for assessing the interactions, binding affinity, and inhibition constants of enterotoxin B when combined with aromadendrin and afzelechin. The structures of aromadendrin (PubChem CID: 122850) and afzelechin (PubChem CID: 442154) were obtained from the PubChem database and were then optimized using OpenBabel software by converting the SDF format to PDB. The crystal structure of enterotoxin B was downloaded from the Protein Data Bank (PDB code: 1SSB, resolution 2.1 Å) and utilized for the MD study. Protein preparation was further executed using AutoDock Tools. All extra residues were isolated from the protein structure, polar hydrogen atoms were included, and water molecules were removed using AutoDock Tools 4.2.2. The energy minimization process was conducted using the GROMACS 2019.6 software, applying the AMBER99SB force field to resolve any steric clashes and to refine the geometry. The structure that underwent energy minimization was then employed for further docking studies. Grid maps were generated with dimensions of 60×60×60 points and a grid spacing of 0.375 Å was chosen across all three axes (X, Y, Z), with the center of this box specified by the coordinates X = 38.565, Y = 50.619, and Z = 44.892. Following this, an autodock run was conducted, featuring 200 docking simulations and 25 million energy evaluations, utilizing the Lamarckian genetic algorithm (LGA). The docking conformation that exhibited the lowest binding energy within the most populated cluster for each system was selected for the molecular dynamics simulation. The interactions between the ligands and enterotoxin B were visualized through the use of VMD1.9.3 and Ligplot+ programs.

Molecular dynamics simulations were conducted to investigate the intricate interactions between enterotoxin B and the ligands aromadendrin and afzelechin. The state-of-the-art GROMACS 2019.6 software, utilizing the robust AMBER99SB force field, was employed to generate detailed topology files for each of the three systems under study. To achieve charge neutrality within the molecular complexes, specific ions, namely chloride (Cl⁻) and sodium (Na⁺), were strategically introduced. The simulations were performed within a defined cubic box fully solubilized with water to mimic a physiological environment. To derive the necessary parameters for all involved ligands, the ACPYPE tool was expertly employed, ensuring a comprehensive representation of the molecular characteristics. Following this preparatory phase, each system underwent a 1-nanosecond (ns) equilibration within both NVT and NPT ensembles, carefully maintaining a temperature of 310 K and a pressure of 1 bar to ensure stability and reliability of the system. Upon completion of the equilibration phase, an extensive MD simulation spanning 200 ns was conducted, utilizing a precise time step of 2 femtoseconds (fs) to capture the dynamic motions of the molecules accurately. The trajectories obtained from these simulations were thoroughly analyzed to provide profound insights into various critical metrics, including the root mean square deviation (RMSD), root mean square fluctuation (RMSF), radius of gyration (RG), solvent accessible surface area (SASA), and the hydrogen bond interactions of enterotoxin B with aromadendrin and afzelechin. This analysis also encompassed the interactions involving the protein and the surrounding solvent, yielding a comprehensive understanding of the molecular dynamics at play.

Figure 1 vividly depicts the distinct binding orientations of aromadendrin and afzelechin within the active site of the enterotoxin B protein, showcasing their intricate molecular interactions. In particular, Figure 1a emphasizes the pivotal van der Waals forces that facilitate the engagement between aromadendrin and enterotoxin B, shedding light on the subtleties of their biochemical relationship. The key amino acids residing in the active site of enterotoxin B were highlighted as TYR198, ASP199, ALA203, PRO202, PHE208, ASP209, TYR213, LYS212, and MET216. This detailed analysis reveals how the hydroxyl group of aromadendrin forms multiple hydrogen bonds with the carbonyl and carboxyl functional groups of ASP 199, the carboxyl group of ASP 209, as well as the carbonyl group of LYS212, illustrating the complex nature of these molecular interactions. Expanding upon these findings, Figure 1b explores the hydrogen bonding dynamics associated with afzelechin and its corresponding active site residues in enterotoxin B. The schematic representation clearly delineates that the hydroxyl groups of afzelechin forge significant hydrogen bonds with the carbonyl and carboxyl groups of ASP 199, the carbonyl group of ASP209, and the carboxyl group of LYS212, further illuminating the sophisticated interplay between these important biochemical partners.

Complementing the visual insights, Table 1 provides an exhaustive overview of the docking parameters, detailing binding energies, inhibition constants, cluster rank, and the number of molecules within each cluster for both aromadendrin and afzelechin as they interact with enterotoxin B. The data indicate that while both compounds exhibit a pronounced affinity for enterotoxin B, aromadendrin stands out with the lowest binding energy, highlighting its superior binding strength in contrast to afzelechin.

Table 1: The obtained docking features predicted by the AutoDock program

Figure 1: The optimal docking pose and molecular interactions of a) Aromadendrin/ enterotoxin B and b) Afzelechin/ enterotoxin B are presented. The figures were generated using VMD version 1.9.3 and Ligplot+ software

Upon the completion of molecular docking computations, we diligently refined the MD simulation calculations to assess the structural stability and residue variations of the complexes. This refinement was primarily conducted to validate the output parameters obtained from the molecular docking studies. The results derived from the MD simulation were utilized to investigate the interactions and structural modifications of enterotoxin B in the presence of aromadendrin and afzelechin. Collectively, these evaluations serve as a foundation for understanding the biological activity in the presence of aromadendrin and Afzelechin, as well as in their unbound forms.

Through RMSD analysis, one can evaluate the conformational changes and stability of complexes throughout the duration of simulations. Figure 2 illustrates the RMSD of free enterotoxin B and its interaction with various ligands. As depicted in Figure 2a, the complex formation of aromadendrin with enterotoxin B results in structural stability when compared to free enterotoxin B. This figure further demonstrates that both the free and bound enterotoxin B systems achieved stabilization after approximately 160 nanoseconds. Conversely, the formation of complexes between afzelechin and enterotoxin B increased structural fluctuations of enterotoxin B, as represented in Figure 2b. According to this plot, the free enterotoxin B attained equilibrium at 160 nanoseconds, while the complex of afzelechin and enterotoxin B equilibrated at approximately 110 nanoseconds. However, as the simulation progresses, the afzelechin complex exhibits greater structural variations in comparison to free enterotoxin B.

Table 2 provides the average values of the molecular dynamics parameters for all four systems over the final 50 nanoseconds. Upon reviewing Figures 2 and Table 2, the binding of aromadendrin and afzelechin significantly impacts the RMSD profile of the enterotoxin B protein. Notably, the average RMSD of enterotoxin B decreased in the presence of aromadendrin, suggesting that aromadendrin may stabilize enterotoxin B. In contrast, the average RMSD of enterotoxin B increased in the presence of afzelechin, which indicates a destabilizing effect on enterotoxin B, as evidenced by a significant alteration in its average RMSD value. Hence, it can be inferred that the binding of afzelechin to enterotoxin B may disrupt its functional capabilities and impede its activity.

Table 2: The average and standard deviations of RMSD, Rg, RMSF, and SASA for free and complexes proteins during the last 50ns

Figure 2: RMSD plots of free and bound systems for a) 1SSB - Aromadendrin, b) 1SSB – Afzelechin

It is essential to undertake a Rg analysis to thoroughly assess the structural compactness and variability of our systems. In our study, we computed the Rg values for both the free protein and its complex with aromadendrin and Afzelechin, enabling us to gain insights into the dynamic structural behavior throughout the simulation period. The findings are visually represented in Figure 3, which clearly illustrates that the binding interactions of aromadendrin and afzelechin lead to a significant reduction in the Rg of enterotoxin B. The mean Rg values, were documented in Table 2, further support our observations. A careful examination of Figure 3 in conjunction with Table 2 reveals that during the latter 50 nanoseconds of the simulation, the average Rg value for enterotoxin B complexed with aromadendrin and afzelechin exhibited a notable decrease. This trend indicates that the association of these ligands profoundly influences the tertiary structure of enterotoxin B, resulting in a marked structural compression. Therefore, we can deduce that the binding of aromadendrin and afzelechin not only promotes structural compression but also induces a degree of destabilization within enterotoxin B, culminating in an overall condensation of the protein's structure.

The RMSF analysis plays a crucial role in assessing the deviations and variability of residues within both free and ligand-bound systems. As illustrated in Figure 4, the RMSF values for the free and complex systems reveal compelling insights. Notably, the majority of amino acid residues in the 1SSB - 1SSB-aromadendrin and 1SSB-afzelechin systems exhibit greater fluctuations compared to free enterotoxin B. This observation suggests that the presence of these ligands significantly influences the dynamic behavior of the protein. Throughout the simulation period, the RMSF values for most residues located in the binding sites tend to increase in the presence of our ligands, indicating enhanced conformational variability.

Furthermore, as presented in Table 2, the mean RMSF value for enterotoxin B in complex with our ligands is notably higher than that of the free enterotoxin B. This finding indicates that when the protein is bound to aromadendrin and Afzelechin, it experiences more pronounced conformational fluctuations in the residual sites compared to its unbound state. Such insights underscore the role of ligand binding in modulating the structural dynamics of the protein, which may have implications for its functional mechanisms.

Figure 3: RG plots of free and bound systems for a) 1SSB - Aromadendrin, b) 1SSB – Afzelechin

Figure 4: RMSF plots of free and bound systems for a) 1SSB - Aromadendrin, b) 1SSB – Afzelechin

The analysis of SASA constitutes a critical methodology for acquiring insights into the surface area of biomolecules that is accessible to solvent, which is essential for understanding molecular interactions throughout the simulation period. Figure 5 illustrates the SASA plot for all experimental systems. The data presented indicate a reduction in the average SASA for the 1SSB - aromadendrin and 1SSB - afzelechin complexes. This decrease can be attributed to the ligands occupying the active cavity of enterotoxin B, thereby interacting with specific amino acids, resulting in a diminished surface area of the protein available to water molecules.

The analyzed systems' average SASA values were computed and summarized in Table 2. This table reveals that the average SASA value for enterotoxin B decreased in the presence of aromadendrin and Afzelechin. Specifically, the SASA value for enterotoxin B complexed with aromadendrin decreased from 133.183 ± 2.303 nm² in its unbound form to 129.525 ± 2.400 nm² in the complex form, while the value for the afzelechin complex with enterotoxin B reduced to 131.266 ± 2.329 nm². Furthermore, the observed decrease in the average SASA for enterotoxin B in the presence of these ligands correlates with a reduction in structural compactness, as indicated by a decrease in the Rg value.

A comprehensive examination of the hydrogen bonds between the protein and its ligands is essential for deciphering the stability of their complexes and the nature of their molecular interactions. As illustrated in Figure 6, throughout a 200 ns simulation period, the dynamics of hydrogen bonding between aromadendrin, afzelechin, and enterotoxin B were analyzed. Our findings reveal that the peak number of hydrogen bonds formed reached five for aromadendrin and four for afzelechin, highlighting the robust stability of these interactions.

Furthermore, Figure 7 demonstrates that the hydrogen bonds established between both ligands and enterotoxin B remained remarkably stable throughout the simulation, reinforcing the notion of complex stability. Notably, the binding of aromadendrin resulted in a notable increase in the average number of hydrogen bonds among atoms within the enterotoxin B protein. This observation suggests that the binding event may induce alterations in the protein's secondary structure, as further substantiated in Figure 7. In contrast, the presence of afzelechin elicited only a minor modification in the hydrogen bonding landscape of enterotoxin B.

Moreover, the interaction between enterotoxin B and solvent molecules exhibited a significant decline in hydrogen bond counts upon the introduction of aromadendrin and afzelechin, as depicted in Figure 8. To provide a quantitative perspective, Table 3 summarizes the average and standard deviations of intra-molecular protein and protein-solvent hydrogen bonds during the final 50 ns of the simulation, thereby corroborating the insights derived from Figures 7 and 8.

To further validate these findings, the Molecular Mechanics Poisson-Boltzmann Surface Area (MMPBSA) method was employed to evaluate the binding energy components, encompassing van der Waals, electrostatic, polar solvation, and SASA energies, utilizing trajectory data from the last 50 ns. Aligning with our docking predictions, the binding free energy calculations obtained via the MMPBSA method unequivocally confirmed the stability of the enterotoxin B complexes with aromadendrin and afzelechin. The consolidated findings, presented in Table 4, illuminate that both complexes establish favorable interactions with enterotoxin B, with aromadendrin demonstrating the strongest binding affinity, thereby corroborating the docking results that indicated the lowest binding energy in the aromadendrin-enterotoxin B complex compared to analogous complexes involving afzelechin.

Figure 5: SASA plots of free and bound systems for a) 1SSB - Aromadendrin, b) 1SSB – Afzelechin

Figure 6: Time dependence of the number of hydrogen bonds for a) 1SSB - Aromadendrin, b) 1SSB – Afzelechin

Table 3: The average and standard deviations of intra-molecular protein and protein-solvent hydrogen bonds during last 50 ns

Figure 7: Protein-Protein H-bond plots of free and bound systems for a) 1SSB - Aromadendrin, b) 1SSB – Afzelechin

Table 4: The average of energies components for complexes analyzed by MMPBSA

Figure 8: Protein-Solvent H-bond plots of free and bound systems for a) 1SSB - Aromadendrin, b) 1SSB – Afzelechin

In the context of existing SEB inhibition research, our molecular docking and MD simulation results for aromadendrin and afzelechin align with and extend prior findings on small molecule interactions with SEB. A study investigated flavonoid binding to SEB, showing that (-)-epigallocatechin-3-gallate (EGCG) binds with higher affinity than other polyphenols and may engage the toxin’s T-cell recognition channel region [41].

Our extensive computational analysis unveils aromadendrin and afzelechin as promising natural inhibitors of enterotoxin B, each exhibiting distinct yet complementary mechanisms of action, aligning to our previous research performed on taxifolin, pachypodol, and camptothecin with enterotoxin B [42]. Our computational binding free-energy analyses similarly position aromadendrin within SEB interaction regions, suggesting potential for functional interference in toxin host engagement. Functional evidence for small molecules comes from studies on tea catechins: a study demonstrated that EGCg neutralizes SEB’s biological activity in vitro and in vivo, dampening cytokine production (TNF-α, IFN-γ) and lethality in mice, thus supporting the translational relevance of polyphenolic SEB binders [43]. Although our results are currently in silico, the MD/MMPBSA stability profiles and energetics position aromadendrin comparably within the class of SEB-binding polyphenols. More potent inhibition strategies have been explored using aptamer antagonists. In Wang et al. study, a SELEX-derived aptamer bound SEB with Kd ≈ 64 nM and effectively blocked SEB-mediated T-cell proliferation and cytokine responses, including in a murine toxic-shock model, indicating functional neutralization at nanomolar affinity [44]. Similarly, the SEB1741 aptamer reduced CD4⁺ T-cell activation and IFN-γ/IL-2 release, demonstrating that nucleic-acid ligands could modulate superantigen effects [45].

These aptamer studies provide experimental benchmarks for affinity and functional inhibition that exceed typical small molecule protein interactions in both specificity and biological activity. At the biologics level, monoclonal antibodies against SEB (e.g., humanized 20B1 variants) exhibit potent toxin neutralization and survival benefit in animal models of intoxication and sepsis, outperforming smaller ligands in affinity and in vivo efficacy [46]. Combinatorial mAb strategies further enhance efficacy by targeting distinct conformational epitopes and altering SEB’s structural dynamics [47]. Relative to these diverse approaches, the natural flavonoids studied here occupy a niche that balances drug-like properties and food safety with moderate binding energetics. While aromadendrin and afzelechin do not yet rival the binding affinities of aptamers or the in vivo neutralization efficacy of monoclonal antibodies, their stable interaction profiles and favorable physicochemical properties make them promising lead compounds for development.

Future work should incorporate biophysical binding assays (e.g., SPR, ITC) and cell-based functional tests (e.g., T-cell activation assays) to quantify affinities and inhibitory potency relative to the referenced benchmarks.

This study provides novel computational evidence that the flavonoids aromadendrin and afzelechin can act as potent inhibitors of SEB, a heat-stable superantigen frequently implicated in staphylococcal food poisoning. By employing molecular docking and MD simulations, we demonstrated that these natural compounds interact strongly with SEB, destabilizing its structure and limiting its toxic potential. From a toxicological perspective, these findings are highly significant, as SEB contamination is a major contributor to foodborne intoxications worldwide, and effective countermeasures remain limited. The ability of dietary flavonoids to inhibit SEB not only highlights their therapeutic promise but also points to potential applications in food safety, nutraceutical development, and toxin mitigation strategies. Future in vitro and in vivo validation of these results will be essential to confirm their protective role and assess their feasibility as natural interventions in preventing foodborne diseases. Overall, this work integrates food science, toxicology, and computational biology to suggest that flavonoids may represent a safe, natural, and sustainable approach to combatting enterotoxin-related food poisoning.

This research was supported by resources supplied by the deputy of financial affairs of the Ghalib University, Kabul, Afghanistan.

The authors declare that they have no competing interests.

The authors would like to express their utmost gratitude to the Board of Directors of Ghalib University, Kabul, Afghanistan for support and motivations, especially Dr. M. I. Noori, Eng. A. Ahadi and Mr. N. A Nadeem. Special thanks to Mr. Mohammad Yousoof Saleh (Head of HR and financial affairs of Ghalib University) and his team, for assisting us in better performance of this study.