Discovery of potential natural therapeutics targeting cell wall biosynthesis in multidrug-resistant Enterococcus faecalis: a computational perspective

Abstract

Background

Identifying therapeutic inhibitors of crucial enzymes involved in the peptidoglycan biosynthesis pathway is pivotal for developing new treatments against multidrug-resistant Enterococcus faecalis V583. MurM, an essential enzyme in this pathway, plays a significant role in the bacterium’s cell wall synthesis, making it an attractive druggable target for novel antimicrobial strategies. This study explored the potential of natural compounds as inhibitors of MurM, aiming to discover promising drug candidates that could serve as the foundation for future therapeutic development.

Methods

The three-dimensional structure of MurM was predicted, optimized, and its binding pocket was analyzed by comparing it with related structures. Over 4,70,000 natural compounds from the COCONUT database were subjected to virtual high-throughput screening (vHTS). The top lead candidates were selected based on their Lipinski’s profile, ADME profile, toxicity profile, estimated binding free energy (ΔG) and estimated inhibition constant (Ki). Interaction pattern analysis was used to evaluate the non-covalent interactions between the inhibitors and key residues in MurM’s binding pocket. Molecular dynamics simulations were performed over 300 ns to assess the structural stability and impact of these inhibitors on MurM’s enzyme.

Results

Three lead compounds—CNP0056520, CNP0126952, and CNP0248480—were identified and prioritized with estimated ΔG ranging from − 9.35 to -7.9 kcal/mol. Molecular dynamics simulations revealed minimal impact on MurM’s overall structure and dynamics, with the candidate inhibitors forming stable protein-ligand complexes. These interactions were supported by several non-covalent interactions between the candidate inhibitors and key residues within MurM’s binding pocket.

Conclusion

These findings suggest that the identified natural product candidates could serve as promising inhibitors of MurM, potentially leading to novel therapeutics targeting cell wall biosynthesis in multidrug-resistant E. faecalis.

Introduction

Enterococci are commonly found in the gastrointestinal tracts of both humans and animals. While most enterococci are non-pathogenic, some are frequently present in artisan fermented foods, which contribute to extending shelf life, enhancing flavor, and improving texture. Notably, enterococci exhibit intrinsic antibiotic resistance and are recognized as opportunistic pathogens in nosocomial or hospital-acquired infections [1,2,3]. They are often isolated from patients with prolonged hospitalization and extensive antibiotic usage. Some Enterococcus strains cause severe infections in humans, including urinary tract infections, surgical wound infections, pneumonia, bacteremia, and bacterial endocarditis [1]. Enterococci commonly show innate resistance to diverse classes of antibiotics [2]. They also acquire antibiotic resistance determinants, including vancomycin resistance, from their environment [3].

E. faecalis and E. faecium are the two most well-known species that cause numerous infections in humans and animals [4]. The bacterium under study is E. faecalis V583, the first vancomycin-resistant clinical isolate reported in the USA [5], whose genome sequencing was already reported [6, 7]. E. faecalis accounts for up to 80% of enterococcal nosocomial infections [8], making it a leading cause of nosocomial infections in Europe and the USA. The increased reports of these bacterial infections are mainly due to the bacterium’s acquired and intrinsic resistance to generally used antibiotics (including vancomycin) and the production of virulence factors involved in adhesion, colonization, resistance, and immune system evasion.

Enterococcal cell walls are primarily composed of peptidoglycan (P.G.), wall teichoic acid (WTA), and lipoteichoic acid (LTA), similar to other Gram-positive bacteria [9]. The P.G. layers form a lattice structure above the phospholipid bilayer, shielding the cell from osmotic stress and pressure. The P.G. structure in E. faecalis features a stem peptide comprising L-alanine, D-glutamine, L-lysine, D-alanine, and D-alanine. Additionally, an L-alanine-L-alanine bridge connects to the ε-nitrogen of the L-lysine side chain in the stem [11, 12]. MurM belongs to the FemXAB family and plays a crucial role in P.G. biosynthesis. MurM is a tRNA-dependent ligase that adds L-alanine to lipid II substrates during cell wall synthesis. This addition is pivotal for properly assembling the peptidoglycan layer, contributing to the structural integrity of the bacterial cell wall [11, 12]. Targeting MurM is a strategic approach in combating antibiotic resistance, as it disrupts cell wall formation and potentially provides selective antimicrobial activity.

Natural compounds have gained significant attention for their diverse and potent biological activities, offering a promising avenue for developing novel therapeutics. The COCONUT database [13] is a pivotal computational resource comprising natural products from a wide range of sources, making it an ideal tool for the comprehensive screening of potential bioactive compounds against various druggable target proteins. In our previous study [14], we explored the structure of MurM and identified potential binding site residues crucial for intermolecular interactions with substrate and ligand molecules. Using the extensive natural compound collection in the COCONUT database, we aim to identify natural candidate inhibitors that target the MurM enzyme. Addressing this research gap could lead to the discovery of novel therapeutics from natural sources targeting cell wall biosynthesis in E. faecalis, providing a new strategy in the fight against antibiotic resistance.

Methodology

Selection of drug target – MurM and structure analysis

The drug target selected for this study is MurM, an enzyme from E. faecalis (strain V583), which plays a crucial role in peptidoglycan biosynthesis. The amino acid sequence of MurM (UniProt Accession Number: Q830W0) [15] was retrieved from the UniProt database. The three-dimensional model of MurM was predicted using a consensus approach and optimized using YASARA [16], as reported previously [14]. Structure analysis was subsequently performed utilizing PyMOL [17] for visualization and PDBsum [18] for detailed structural insights concerning fold, architecture, topology, secondary structural class and secondary structural compositions. The possible binding pocket of MurM was analyzed by comparing its related structures. We have performed comparative protein sequence alignment and structural superposition of MurM and its related proteins using the Multalin (https://multalin.toulouse.inra.fr/multalin/) and PyMOL (https://www.pymol.org/) programs to identify and prioritize the binding site residues.

In house library preparation from the COCONUT database and virtual high throughput screening

To identify potential natural therapeutics targeting MurM, we utilized the COCONUT database [13], which comprises natural products from diverse sources. The molecules downloaded from COCONUT were initially filtered using the DataWarrior [19] program by applying Lipinski’s rule of five. This rule assesses drug-likeness based on various physiochemical parameters, including hydrogen bond donors (HBD) and hydrogen bond acceptors (HBA), lipophilicity (LogP), molecular weight (M.W.), and the number of rotational bonds. Additionally, toxicological properties, including carcinogenicity, mutagenicity, reproductive toxicity, and irritation potential, were evaluated. The filtered compounds were further analyzed using SwissADME [20] to predict absorption, distribution, metabolism, and excretion (ADME) profiles. Parameters such as aqueous solubility (AS), blood-brain barrier penetration (BBB), CYP2D6 inhibition, topological polar surface area (TPSA), and human intestinal absorption (HIA) were considered. The promising lead or drug-like candidates from natural sources were prioritized based on their binding orientation within the MurM’s active site, key interacting residues, and the number of intermolecular hydrogen bonds (H-bonds) or other non-covalent interactions formed with active site residues. To identify and prioritize potential lead candidates against the MurM protein from E. faecalis, we conducted virtual high-throughput screening (vHTS) in triplicates using AutoDock Vina [21] against an in-house library of natural compounds from the COCONUT database. The docking parameters for the virtual high-throughput screening (vHTS) were configured as follows: the center coordinates were set at X = 8.6 Å, Y = 29.3 Å, and Z = 20.6 Å. The grid box was defined based on the binding site residues of MurM (Trp31, Lys35, Trp38, Arg215, Tyr219, and Tyr315), with dimensions of 30.9 Å along the X-axis, 26 Å along the Y-axis, and 34 Å along the Z-axis. The exhaustiveness parameter was set to 10, and the energy range was adjusted to 4.

Molecular docking studies of prioritized lead candidates with MurM

We have additionally conducted molecular docking calculations of MurM against specific prioritized lead candidates to cross validate the vHTS results. We used molecular graphics laboratory (MGL) tools to prepare protein and ligand molecules. The ligands were in the PDBQT format, which contains XYZ coordinates (PDB), partial charges (Q), and atom types (T). One of the ligands, called Control, as Deoxyribofuranosephosphate, was chosen as an inhibitor based on prior study [8]. The AutoDock [22] was used to conduct molecular docking calculations by leveraging the possible binding sites of MurM. The molecular docking calculations (site-specific) enable an understanding of the binding orientation, binding mode, ∆G, Ki, and various intermolecular interactions between protein and ligand molecules. The grid box was defined based on the binding site residues of MurM (Trp31, Lys35, Trp38, Arg215, Tyr219, and Tyr315), with dimensions (npts) set to X = 74, Y = 68, and Z = 62. The center coordinates were X = 10.86 Å, Y = 27.53 Å, and Z = -20.60 Å. The Lamarckian Genetic Algorithm (LGA) was employed for the MurM-ligand docking, with each ligand undergoing 1000 docking runs and a population size of 300. The maximum number of generations and energy evaluations were 27,000 and 2,500,000, respectively, while the mutation and crossover rates were maintained at 0.02 and 0.8. After analyzing the MurM-ligand complexes using MGLTools 1.5.7 and PyMOL [17] programs, we selected those with favorable estimated free energy of binding values for further molecular dynamics (MD) simulations. BIOVIA Discovery Studio was also utilized to analyze and visualize the protein-ligand complexes, enabling us to investigate the intermolecular interactions, interaction patterns, interacting residues, and types of intermolecular interactions in detail.

Molecular dynamics simulations

The MD simulations provided critical insights into the structural stability and dynamic behavior of the MurM-ligand complexes and the unbound form of MurM, facilitating the identification of potential lead molecules for further experimental studies. These simulations were conducted using GROMACS (GROningenMAchine for Chemical Simulations) Version 2021.3 [23]. Initially, the input files in GROMACS format (GRO), along with the topology (TOP) and positional restraint (POSRE) files, were generated from the MurM’s PDB file. A cubic box was created around both the unbound MurM and its ligand-bound complexes (System Size: 9.241 nm, 6.042 nm, 5.882, Diameter: 9.957 nm, Center: 2.182 nm, 2.357 nm, -2.351 nm, Distance (distance between centre of the protein and edge of the box): 1.0 nm). The TIP4P water model was employed to solvate the systems by adding water molecules to the simulation box. In the present study, we chose the CHARMM27 [31] force field due to its well-established accuracy in modeling protein structural dynamics, making it particularly suitable for the MurM protein from E. faecalis. CHARMM27 is specifically optimized for protein-ligand interactions, which is essential for studying binding mechanisms in MurM. The TIP4P water model [32] was selected for its ability to accurately reproduce water’s physical properties, such as hydrogen bonding and dielectric behavior, which are critical for simulating realistic solvation effects. Thus, the CHARMM27 force field combined with the TIP4P water model provides a well-suited and validated framework to explore the structural and dynamic properties of the MurM protein from E. faecalis, as well as its interactions with ligands. Counter ions were added when necessary to ensure charge neutrality, except when the systems were already neutral. Since GROMACS could not generate topology files for the ligand molecules, so we used the SwissParam web server [24] to construct the selected ligands’ topology (ITP) and provide PDB files. The gmx editconf tool was then used to convert the ligand PDB files to GRO files. Subsequently, the MurM systems (unbound and ligand-bound complexes) underwent energy minimization using the steepest descent algorithm to stabilize them and remove the larger strains or other geometrical errors. The systems were then equilibrated in two phases: the NVT equilibration (constant volume and temperature) and the NPT equilibration (constant pressure and temperature), each for 3 ns with aid of POSRE files, which are generated during the gmx pdb2gmx step. A velocity-rescaling thermostat and a Parrinello-Rahman barostat were used to maintain the desired temperature and pressure. Following equilibration, production MD simulations were carried out for 300 ns for all MurM systems. During the production MD steps, we removed the positional restraints to observe the natural dynamics of the proteins and ligand-bound complexes. The molecular dynamics parameters (MDP) files were adjusted for the ligand-bound MurM systems during energy minimization, NVT equilibration, NPT equilibration, and production MD steps. The resulting MD simulation trajectories were corrected by removing periodic boundary conditions (PBC). These corrected trajectories of the unbound and ligand-bound MurM complexes were then analyzed using various pivotal structural parameters, including RMSD, RMSF, Rg, SASA, intermolecular hydrogen bonds, and Principal Component Analysis (PCA). Specific GROMACS commands such as gmx rms, gmx rmsf, gmx gyrate, gmx sas, gmx hbonds, gmx covar, and gmx anaeig were used for comparative structural analysis of MurM and ligand-bound MurM systems. In addition to stability analysis, we conducted MM/PBSA-based binding free energy and individual residue energy decomposition analysis for ligand-bound MurM complexes. We used the converged trajectories ranging from 1 to 150 ns, with snapshots taken at 250 ps intervals. The g_mmpbsa tool [33] calculated the binding free energy values of the selected ligand-bound MurM complexes. The binding free energies were calculated using the following equation:

∆G_binding = ∆_vdw + ∆_ele + ∆_pol + ∆_np.

Where ∆_vdw, ∆_ele, ∆_pol and ∆_np represent the contributions of van der Waals energy, electrostatic energy, polar solvation energy, and non-polar solvation energy to the binding free energy (∆G_binding).

Results

MurM structure and its binding pocket

MurM from E. faecalis does not have its three-dimensional structure in the Protein Data Bank. Therefore, the structure of MurM (Fig. 1) was predicted using a consensus approach with the aid of multiple in silico tools and was optimized and reported previously [14]. The structure validation results of MurM were presented in Supplementary Data. The MurM model is structured as a compact unit with a globular domain comprising two sub-domains. The binding site for the substrate is situated at the interface between these two sub-domains of MurM. Each sub-domain is characterized by two twisted β-sheet cores surrounded by α-helices. Subdomain 1A encompasses amino acid residues 1–153 and 382–401, while sub-domain 1B consists of amino acid residues 154–241 and 294–381. In addition, the coiled-coil helical arm spans residues 242–293. A central mixed polarity sheet comprises five strands within each sub-domain, accompanied by four α-helices. Specifically, two α-helices are positioned above the sheet, parallel to the β strands. Overall, as per the results of the PDBSum [18] generate tool, the 3D structure of MurM contains 16 strands, 18 α-helices, 18 β turns, and two βαβ units. A core mixed polarity sheet comprising five strands and four α-helices is present in every sub-domain.

Fig. 1

Three-dimensional structure of MurM with colors representing different domains of the MurM protein

Using the Position Specific Iterative - Basic Local Alignment Search Tool (PSI-BLAST) [25] to compare MurM against the PDB database, we identified the structures of seven closely related proteins: (a) 6SNR, (b) 1XE4, (c) 1XF8, (d) 1LRZ, (e) 4II9, (f) 1NE9 and (g) 3GKR. Subsequently, comparative sequence alignment and structural superposition were carried out on MurM and its related proteins. Both MurM and its homologs exhibit similarities in substrate-binding site residues. The predicted model of MurM highlights potential binding site residues, including Trp31, Lys35, Trp38, Arg215, Tyr219, and Tyr315, which are highly conserved for MurM and its related proteins (Fig. 2). The comprehensive results of the comparative sequence alignment and structural superposition are provided in the Supplementary Data.

Fig. 2

Three-dimensional structure of MurM highlighting the potential active site residues identified through structural comparison technique with related proteins

Virtual high throughput screening, physicochemical properties and ADMET analysis

To identify possible natural candidate inhibitors against MurM, a vHTS was performed on the COCONUT Database, which comprised 4,72,270 natural molecules. Molecular docking attempts to anticipate the optimal binding orientation of candidate inhibitors into the binding pocket of the druggable target proteins. In the current study, we employed the structure-based virtual screening (SBVS) technique, a vital tool for elucidating protein-ligand interactions and identifying potential drug-like candidates. To inhibit the function of the target protein - MurM, the molecules obtained using vHTS should be oriented in the aforementioned binding pocket (interface between two domains). We chose Deoxyribosefuranosephoaphate as the control molecule, based on the previous study [8], with an IC₅₀ value of 100 µM. The systematic workflow for identifying promising candidate inhibitors is shown in Fig. 3.

Fig. 3

Workflow outlining the systematic process to identify and prioritize natural lead candidates against MurM using the COCONUT database

After retrieving 4,72,270 molecules from the COCONUT database, we further obtained 98,782 lead molecules with significant drug-like properties. These ligands were assessed based on criteria where drug-likeness (M.W. ≤ 500 Da, HBD ≤ 5, HBA ≤ 10 and LogP ≤ 5), rotatable bonds (≤ 10), and toxicity (none) were required and among the initial set comprising 98,782 molecules, only 3,285 molecules exhibited favorable ADME profiles, significant pharmacokinetic properties, and no detectable side or adverse effects. For SwissADME pharmacokinetic analysis, we selected ligands with the following characteristics: HIA + and BBB+ (able to penetrate the blood-brain barrier), solubility in aqueous medium, TPSA value less than 100 Ų, and non-inhibition of CYP2D6 (drug-metabolizing enzyme) isoforms. Subsequently, the top three ligand molecules were chosen based on their orientation in the MurM binding site. From the vHTS technique, we selected three ligands (CNP0056520, CNP0126952, and CNP0248480) positioned within the MurM binding site for further molecular docking and MD simulations studies. These compounds exhibit estimated free energy of binding values (ΔG) obtained from AutoDock of -8.5 kcal/mol, -7.93 kcal/mol, and-9.35 kcal/mol, whereas the control molecule exhibited − 6.44 kcal/mol. Tables 1, 2, 3 and 4 provide details on the estimated binding affinity, drug-likeness, ADME, and toxicity profiles of prioritized lead candidates against MurM.

Table 1 Prioritized lead candidates and control molecule with ΔG, Ki, molecular formula, and chemical structure

Table 2 Toxicity profiles of prioritized lead candidates along with control molecule

Table 3 ADME profiles of prioritized lead candidates along with control molecule (HIA Human Intestinal Absorption, BBB blood-brain barrier, TPSA Total Polar Surface Area and CYP2D6 Cytochrome)

Table 4 Lipinski’s profiles of prioritized lead candidates along with control molecule (M.W. molecular weight, clogP Consensus log P octanol-water partition coefficient, LogP lipophilicity, HBA Hydrogen Bond Acceptor, HBD HydrogenBond Donor

Interaction studies and interaction pattern analysis of MurM-ligand complexes

The three prioritized lead candidates and control molecule were docked against the MurM protein using AutoDock software. The docking orientation of the proposed lead candidates along with control molecules is illustrated in Fig. 4. CNP0056520 exhibited an estimated inhibition constant (Ki) of 590.43 nM and an estimated free energy of binding (ΔG) of -8.5 kcal/mol with MurM. As illustrated in Figs. 5 and 6, analysis of MurM-ligand interactions revealed that the ligand formed hydrogen bonds, pi interactions, and Van der Waals interactions with key residues within the MurM binding site (Table 5). Specifically, the ligand displayed hydrogen bond interactions with Ile313 (2.39 Å), Leu314 (2.85 Å), pi-alkyl interaction with Trp31 (4.77 Å), and Ala316 (5.37 Å). Furthermore, the molecule exhibited Van der Waals interactions with Leu26, Lys35, Gln153, Pro154, Tyr219, Glu312, Tyr315, Trp332, Met349, Gly350, Gly351, and Leu360, contributing to the stabilization of the MurM-ligand complex. To visually represent MurM-ligand interactions, two-dimensional and three-dimensional ligand interaction diagrams were generated using Discovery Studio and PyMOL tools, respectively, providing a graphical illustration of MurM’s spatial organization and interactions with the CNP0056520, enhancing understanding of the interaction pattern and its binding mechanism towards the binding pocket of MurM.

Fig. 4

Molecular docking orientation of prioritized lead candidates and control molecule towards the substrate binding site of MurM

Fig. 5

Electrostatic surface potential map analysis of prioritized lead candidates and control molecule towards substrate binding cleft of MurM

Fig. 6

Two-dimensional schematic ligand interactions depicting prioritized lead candidates and control molecule binding to MurM

The MurM_CNP0126952 complex showed an estimated ΔG of -7.93 kcal/mol and a Ki value of 1.53 µM. The ligand formed specific interactions with distinct residues in this complex, detailed in Table 5. Hydrogen bonding interactions were observed with the Leu314 (2.73 Å) residue of MurM. Additionally, pi-alkyl interactions with Leu26 (4.01 Å) and Tyr315 (4.59 Å, 5.30 Å) were observed. Furthermore, Van der Waals interactions with Trp31, Lys35, Leu214, Tyr219, Ile313, Ala316, Trp332, Met349, Gly350, Gly351, Lys360, and Lys364 amino acid residues contributed to the binding stability of the complex.

The MurM_CNP0248480 complex exhibited a calculated ΔG of -9.35 kcal/mol, with an estimated Ki value of 140.06 nM. Hydrogen bonding interactions were observed between the ligand and Gly351(2.98 Å) residue, crucial for the MurM-ligand interaction. The ligand’s pi-alkyl interaction with Leu214, Tyr315, Ala316 and Trp332 residues further stabilized the MurM_CNP0248480 complex. Additionally, Van der Waals contacts were identified with Trp31, Met148, Gln153, Pro154, Ile212, Ser213, Glu312, Ile313, Leu314, Met349, Glu350, Leu360, and Lys364 residues, enhancing the overall stability of the MurM complex.

In contrast, the MurM-Control complex exhibited an estimated ΔG value of -6.44 kcal/mol and a Ki value of 19.03 µM. This compound interacted with several key amino acid residues, including Lys35, Trp38, Arg215, and Tyr219, catalytically conserved in MurM and present in related structures. Hydrogen bonding interactions were observed between the ligand and MurM residues Asn37 and Lys60 (Table 5), significantly influencing the ligand’s binding affinity to the target protein. Additionally, pi interactions were detected with MurM’s Lys35 residue, contributing to the stabilization of the MurM-Control complex. Furthermore, Van der Waals interactions with Pro61, Arg215, and Tyr219 amino acid residues further enhanced the overall stability of the complex.

The investigation of MurM binding with the prioritized natural candidate inhibitors revealed that all three compounds engaged in intermolecular interactions with binding site residues, exhibiting ΔG values below − 6.44 kcal/mol. Moreover, aromatic rings increased hydrophobic interactions for these ligands, interacting with MurM’s binding site residues such as Trp31, Lys35, Tyr219, and Gly351. These interactions, including pi and Van der Waals contacts with conserved amino acid residues, underscored their strong binding affinity with MurM. Such intermolecular interactions are crucial in identifying promising lead candidates for the drug design and development process. Additionally, electrostatic surface potential map analysis (Fig. 6) demonstrated that MurM’s substrate-binding pocket predominantly comprises positively charged amino acid residues. Furthermore, all four ligands, including the Control, were appropriately oriented within the MurM substrate-binding cleft, as depicted in Fig. 6.

Table 5 Various intermolecular interactions of prioritized lead candidates and control molecule towards MurM

Molecular dynamics simulations of MurM and its complexes

The molecular dynamics (MD) simulations were performed to evaluate the structural impact of several ligand molecules on MurM protein, such as stability, flexibility, integrity, compactness, and folding behavior. These simulations used the unbound MurM protein and MurM-ligand complexes, including a control molecule and three potential and prioritized natural candidate inhibitors: CNP0056520, CNP0126952, and CNP0248480. The results of these MD simulations are illustrated in Figs. 7, 8 and 9 and average values obtained from MD simulations in Table 6 providing a comprehensive view of the protein’s behavior under the influence of various lead molecules.

Fig. 7

Root mean square deviation (RMSD), Radius of Gyration (Rg), Hydrogen Bonds, and Solvent Accessible Surface Area (SASA) analyses for MurM and its ligand-bound complexes (Color codes: Black–MurM, Red–Control, Green–MurM_CNP0056520 complex, Blue–MurM_CNP0126952, Yellow–MurM_CNP0248480)

Fig. 8

Root mean square fluctuation (RMSF) analysis of MurM and its complexes (Color codes: Black–MurM, Red–Control, Green–MurM_CNP0056520 complex, Blue–MurM_CNP0126952, Yellow–MurM_CNP0248480)

Fig. 9

Essential dynamics analysis of MurM and its complexes (Color codes: Black–MurM, Red-Control, Green- MurM_CNP0056520 complex, Blue-MurM_CNP0126952, Yellow-MurM_ CNP0248480)

The RMSD analysis, shown in the top-left panel of Fig. 7, was used to assess the structural stability of MurM and its complexes throughout the 300 ns simulation period. The RMSD plot reveals that the unbound MurM protein, represented by the black line, stabilized quickly and maintained a relatively low average RMSD value of 0.34 nm. This indicates that the native MurM structure is inherently stable under the MD simulation conditions. In contrast, the MurM-ligand complexes exhibited slightly higher average RMSD values, suggesting that ligand binding induces some degree of structural perturbation or changes. Among the complexes, MurM_CNP0126952 (blue line) showed the highest average RMSD of 0.49 nm, indicating that this ligand caused the most significant structural changes in MurM. The MurM_CNP0248480 complex (yellow line) followed closely with an average RMSD of 0.48 nm. The MurM_CNP0056520 complex (green line) and the control complex (red line) showed intermediate RMSD values of 0.46 nm and 0.40 nm, respectively. These results suggest that while ligand binding significantly affects MurM’s structure, the changes are relatively modest, with all average RMSD values remaining below 0.5 nm. This implies that the overall fold of the protein is primarily maintained upon ligand binding despite some local structural adjustments.

The RMSF analysis, presented in Fig. 8, provides insights into the flexibility of different regions of MurM protein and its complexes. The RMSF plot reveals that the C-terminus of MurM, approximately residues 350–400, displays the highest degree of flexibility across all systems, with fluctuations reaching 1–1.2 nm. This high flexibility in the C-terminal region could be significant for the MurM function, possibly allowing for conformational or structural changes necessary for the substrate binding or catalytic activity. Additionally, notable fluctuations were observed in the amino acid residues 100–150 and 250–275 of the MurM systems. These regions likely correspond to loop structures in the MurM protein, which often exhibit higher flexibility. The increased flexibility in these regions could be vital for the protein’s dynamics and potentially for its interaction with substrate or other ligand molecules. Generally, the RMSF patterns of the MurM-ligand complexes showed varying degrees of fluctuation. Among all the systems, the MurM_CNP0126952 complex exhibited the highest average RMSF of 0.24 nm, indicating that this ligand increases the overall flexibility of MurM. In contrast, the unbound MurM and the MurM_CNP0056520 complex showed the lowest average RMSF of 0.18 nm, suggesting these two systems are the least flexible. The control and MurM_CNP0248480 complexes showed intermediate flexibility with average RMSF values of 0.21 nm and 0.23 nm, respectively. These differences in RMSF patterns indicate that ligand binding can modulate the flexibility of specific regions in MurM, which could affect the protein’s function and its interactions with substrate and ligand molecules.

The Rg analysis, illustrated in the top-right panel of Fig. 7, provides information about the overall compactness, folding properties and shape of the MurM protein and its complexes. The Rg plot shows that the binding of the chosen ligands had minimal influence on the compactness of the MurM structure. The average Rg values for all MurM systems fell within a narrow range of 2.65–2.73 nm, indicating that the overall shape and packing of the protein remained largely unchanged upon the influence of ligand binding. Specifically, the control complex showed the lowest average Rg of 2.65 nm, suggesting it might induce a slightly more compact structure. On the other hand, the MurM_CNP0248480 complex exhibited the highest average Rg of 2.73 nm, indicating a marginally looser packing. The unbound MurM, MurM_CNP0056520 complex, and MurM_CNP0126952 complex showed intermediate Rg values of 2.67 nm, 2.68 nm, and 2.71 nm, respectively. These slight variations in Rg values suggest that while ligand binding may cause subtle changes in the protein’s compactness, it does not significantly alter the overall shape and folding behavior of MurM.

The SASA analysis, presented in the bottom-right panel of Fig. 7, was used to determine the surface area of MurM exposed to the solvent in its free form and ligand-bound complexes. The results confirmed that ligand binding did not significantly alter the overall structure or solvent exposure of the MurM complexes. The average SASA values ranged from 224.0 nm² to 229.629 nm², indicating only minor changes in solvent exposure upon the influence of ligand binding. The unbound MurM and the control complex showed the lowest SASA values of 224.6 nm² and 224.0 nm², respectively, suggesting these structures have the least exposed surface area. The MurM_CNP0126952 complex exhibited the highest SASA of 229.629 nm², indicating that this ligand might induce a slight expansion of the protein surface. The MurM_CNP0056520 and MurM_CNP0248480 complexes showed intermediate average SASA values of 228.89 nm² and 228.012 nm², respectively. These results suggest that ligand binding causes slight changes in MurM’s surface exposure, potentially affecting its interactions with the solvent and other molecules in its environment.

The analysis of intermolecular hydrogen bonds, shown in the bottom-left panel of Fig. 7, provided valuable information about the nature and strength of the interactions between the ligand molecules and MurM protein. Interestingly, the MurM_CNP0126952 complex and the control molecule formed the most stable hydrogen bonds, with an average of 2 hydrogen bonds each throughout the MD simulation. This suggests that these ligands form specific and relatively strong interactions with the MurM protein via hydrogen bonding interactions. In contrast, the MurM_CNP0056520 complex formed the fewest hydrogen bonds, with an average of 0.42 (rounded to 0), indicating weaker hydrogen bonding interactions. The MurM_CNP0248480 complex showed an intermediate level of hydrogen bonding, with an average of 1.25 (rounded to 1) hydrogen bonds. These differences in hydrogen bonding patterns suggest varying modes of intermolecular interaction between MurM and the different ligands, which could affect their binding affinity and potential inhibitory activity.

The Principal Component Analysis (PCA), illustrated in Fig. 9, provides insights into the essential motions of the MurM systems. The 2D projection of the trajectory plot shows distinct cluster distributions for each system, represented by different colors. Notably, the ligand-bound complexes, especially those represented by yellow, blue, and red clusters, occupy larger areas in the projection of conformational space than the unbound MurM (black cluster). This indicates that ligand binding increases the conformational and essential dynamics of MurM.This observation is further supported by the trace of covariance matrix values presented in Table 6. The ligand-bound complexes show higher trace values (361.351 to 584.863 nm²) than the unbound MurM (345.85 nm²). Among the complexes, MurM_CNP0126952 exhibits the highest trace value of 584.863 nm², indicating that this ligand induces the most significant changes in the essential motion of the MurM structure. These results suggest that ligand binding affects local structural features and influences the essential dynamics of MurM protein, potentially altering its functional behavior.

Table 6 Time-averaged structural properties obtained from MD simulations of MurM and its complexes

In addition to global and essential dynamics analyses, we conducted MM/PBSA-based binding free energy calculations and individual residue energy decomposition analysis on the ligand-bound MurM complexes. This allowed us to gain deeper insights into the binding affinity behavior and to identify the key residues responsible for interactions with the MurM protein. These results provide a deeper understanding of the binding affinity between MurM and the proposed lead compounds, offering valuable insights into their potential as inhibitors. The control complex, with a positive binding energy of 45.34 kcal/mol, indicates a weak or unfavorable interaction with MurM, suggesting that it may not be an ideal candidate for binding under the conditions analyzed. Among the tested ligands, CNP0126952 exhibited the most favorable binding energy of -19.29 kcal/mol, indicating a strong and energetically favorable interaction with the MurM protein. This significant negative binding energy highlights its potential as a potent ligand with promising inhibitory effects. CNP0056520, with a binding energy of -17.24 kcal/mol, also demonstrated a strong interaction, although slightly weaker than CNP0126952. This suggests it could still be a viable candidate, offering a reasonable binding affinity to MurM. CNP0248480, with a binding energy of -14.77 kcal/mol, displayed the weakest interaction among the three ligands. While its binding energy is still negative, indicating a favorable interaction overall, its lower value suggests that its binding affinity to MurM is comparatively less potent than the other two ligands. The ranking of these ligands based on their binding affinities is as follows: CNP0126952 > CNP0056520 > CNP0248480. These quantitative insights help prioritize these lead compounds and provide a rational basis for further experimental validation. Such findings could drive the next steps in drug development efforts to inhibit MurM, a critical protein in bacterial cell wall biosynthesis.

The individual residue energy decomposition analysis for the active site residues of MurM across three complexes indicates varying binding affinity effects. Trp31 shows consistently stabilizing contributions with values of -0.4358, -0.8003, and − 0.2076 for all three complexes, suggesting it plays a positive role in ligand interaction. In contrast, Lys35 exhibits highly destabilizing effects, particularly in the MurM_CNP0056520 complex with a value of 13.1327, alongside mildly destabilizing values in the MurM_CNP0126952 complex and MurM_CNP0248480 complex (0.6016 and 0.3271), indicating it may inhibit binding. Residue Trp38 has minimal impact, with values close to neutral across all complexes (0.1283, -0.0236, and 0.0062), suggesting it does not significantly influence ligand interactions. Residue Arg215 shows a notable destabilizing effect in the MurM_CNP0056520 complex (12.435) but stabilizing contributions in MurM_CNP0126952 and MurM_CNP0248480 (-0.2140 and 0.1817) complexes, implying its role varies significantly depending on the environment. Similarly, Tyr219 presents a mix of effects, being mildly destabilizing in MurM_CNP0056520 (0.2058) but stabilizing in MurM_CNP0126952 (-0.3253), while remaining neutral in MurM_CNP0248480 complex (-0.0638). Finally, Tyr315 is predominantly stabilizing, especially in MurM_CNP0126952 with a significant value of -1.9576, while it shows mild destabilizing contributions in MurM_CNP0056520 and MurM_CNP0248480 (0.2360 and − 0.7112) complexes. The detailed results of individual residue energy decomposition analysis were given in supplementary data.

Table 7 MM/PBSA-based binding free energy results of various MurM-ligand complexes (only the three prsomising lead candidates that exhibited negative binding energy values were presented here)

Discussion

The comparative analysis of the three natural candidate inhibitors (CNP0056520, CNP0126952, and CNP0248480) with the control molecule [8] reveals interesting structural insights into their potential as prioritized lead candidates. All three compounds demonstrated significantly improved binding affinities compared to the control, as indicated by their more negative estimated ΔG values and lower Ki values. This enhanced binding affinity is a promising sign of their potential effectiveness as candidate inhibitors.

CNP0126952 has the strongest binding free energy among the three candidates, exhibiting a binding energy of -19.29 kcal/mol (Table 7). This compound also displays the highest flexibility, with elevated average RMSD, RMSF, and covariance matrix trace values, suggesting its adaptability to different binding modes and conformations. Moreover, it forms two hydrogen bonds with the MurM protein, which could contribute to its binding stability.

CNP0056520 also has significant binding free energy at -17.24 kcal/mol, demonstrating the most stable behaviour in molecular dynamics attributes. Its lower average RMSF and covariance matrix trace values indicate a more rigid interaction with the MurM protein. It could be advantageous for maintaining a consistent binding pose but may limit adaptability to conformational changes in the protein structure. This compound forms one hydrogen bond with the target MurM protein, contributing to its structural stability.

CNP0248480 has a binding free energy of -14.77 kcal/mol, still presents a promising combination of strong binding affinity and balanced molecular dynamics properties. The moderate increases in RMSD, RMSF, and Rg values suggest that this compound may undergo slight conformational changes in the MurM structure upon binding, which could enhance its molecular mechanism. It also forms one hydrogen bond with MurM, which may contribute to its binding stability.

In addition to performing global and essential dynamics analyses, MM/PBSA-based binding free energy calculations and individual residue energy decomposition analyses were conducted on the ligand-bound MurM complexes. These analyses helped identify the key residues involved in binding interactions and provided insights into the binding affinities of various ligands. Among the ligands, CNP0126952 showed the most favorable binding free energy, suggesting its strong inhibitory potential, while CNP0056520 and CNP0248480 also demonstrated promising interactions, although with slightly lower binding affinities. The individual residue energy decomposition analysis highlighted key residues like Trp31 and Tyr315 as stabilizing factors, while Lys35 and Arg215 were found to have variable, often destabilizing effects.

When considering the overall profile of prioritized natural candidate inhibitors, CNP0248480 appears to offer the most promising combination of good binding affinity and balanced molecular dynamics attributes. Its ability to form a hydrogen bond while maintaining flexibility suggests that it could effectively bind to the MurM protein while also accommodating potential molecular mechanisms. However, the selection of the best natural candidate inhibitors may ultimately depend on the specific requirements of drug design and discovery. If the target protein requires a more flexible ligand for optimal interaction, CNP0126952 might be the ideal choice. Conversely, if protein-ligand complex stability is a priority, CNP0056520 could be more suitable.

Therefore, the combined results from virtual high-throughput screening, molecular docking, and MD simulations suggest that these three natural candidate inhibitors are exceptionally well-suited for drug design and discovery targeting the MurM protein from E. faecalis. This organism is responsible for nosocomial infections and demonstrates resistance to diverse classes of antibiotics, making these findings particularly significant in the fight against antibiotic resistance [26,27,28,29,30]. Previous studies have effectively utilized structure-based virtual screening and molecular dynamics approaches to identify inhibitors for essential proteins in Staphylococcus aureus [34,35,36,37], providing a solid foundation for our investigation into the MurM protein. It is important to note that while these computational hypotheses provide valuable structural insights, they should be validated through further experimental methods. Further investigations, such as in vitro binding assays, functional assays, and X-ray crystallography studies, would be crucial to confirm the predicted interactions and effects on the target MurM protein. These additional studies would provide a more comprehensive understanding of how these three natural candidate inhibitors interact with the target MurM protein in a biological context, ultimately guiding the selection of the most promising candidates for further clinical trials.

Conclusion

The present computational study successfully identified three natural lead candidates as potential inhibitors of the MurM enzyme in E. faecalis through virtual high throughput screening, molecular docking, and MD simulations. The lead compounds demonstrated strong binding affinities and stable interactions with catalytically important residues in MurM’s structure. Notably, the binding of these compounds did not significantly alter the overall structure, flexibility, or dynamics of MurM, suggesting they may hinder enzyme function without disrupting its structural integrity. The consistent performance of these three compounds across multiple analyses, including RMSD, RMSF, Rg, SASA, and hydrogen bonding, further supports their promise as lead candidates. While experimental validation is needed, this study provides a strong foundation for developing potential molecules targeting cell wall synthesis in multidrug-resistant E. faecalis. Using natural compounds may offer advantages in terms of reduced toxicity and novel mechanisms of action. Future work should optimise these lead compounds and evaluate their binding efficacy against resistant E. faecalis strains in vitro and in vivo.