Molecular Modeling and simulation-based identification of inhibitors against new Delhi Metallo-Lactamase 1: Implications for bacterial antibiotic resistance

1 Introduction

Emergence of bacterial antibiotic resistance has become a major clinical health concern (Ahmed et al., 2021; Faheem et al., 2013; Muteeb et al., 2017). β-lactamases enzyme is one of the more recently discovered agent in the bacterial arsenal for antibiotic resistance, and can efficiently hydrolyze β-lactam antibiotics such as penicillins, cephalosporins, monobactams, and carbapenems (Wang et al., 2021). β-lactamases function by hydrolyzing the β-lactam ring available in β-lactam antibiotics. The hydrolysis reaction involves the nucleophilic attack on the carbonyl carbon of the β-lactam ring, leading to the opening of the ring and loss of antibiotic activity. The production of β-lactamases by bacteria significantly complicates the treatment of bacterial infections. Infections caused by organisms that produce β-lactamase are linked to higher morbidity and mortality rates, extended hospital stays, and elevated healthcare costs. The presence of β-lactamases, particularly extended-spectrum beta-lactamases (ESBLs) and carbapenemases like NDM-1, in pathogenic bacteria such as Escherichia coli and Klebsiella pneumoniae has led to outbreaks in healthcare settings and community environments.

Protein homology and underlying molecular mechanisms of action can be used to classify β-lactamases, as suggested by Ambler, into 4 different groups (A, B, C & D). Class A β-lactamases are “serine β-lactamases”, which include enzymes like TEM and SHV, often found in Gram-negative bacteria. Class C β-lactamases are also serine-based, commonly referred to as AmpC enzymes, which are chromosomally encoded and inducible. Class D β-lactamases, known as oxacillinases, have a serine residue at their active site and are capable of hydrolyzing oxacillin and other β-lactams. NDM-1 belongs to class B, also known as metallo-β-lactamases (MBLs). Unlike other classes, metallo-β-lactamases (MBLs) need divalent metal ions, usually zinc, for their enzymatic activity and can break down a broad spectrum of β-lactam antibiotics, including carbapenems and cephalosporins (Hall and Barlow, 2005). Furthermore, Bush and Jacoby have classified β-lactamases into three functional groups (1, 2 and 3). Cephalosporinases are in group 1, while penicillinases, oxacillinases, serine-based carbapenemases, and ESBLs are in group 2. Metal ion-requiring carbapenemases are classed into group 3 (Bush and Jacoby, 2010).

NDM-1, belonging to Class B1 (Ambler’s classification) and Group 3a (Bush and Jacoby’s classification), is a very potent β-lactamase which can inactivate β-lactam antibiotics, including cephalosporins, carbapenems and penicillins (Brem et al., 2016)7. Mechanistically, NDM-1 selectively hydrolyzes amide bonds in β-lactam antibiotics. NDM-1 was discovered in a Swedish patient in India in 2008, who was diagnosed with urinary tract infection of Klebsiella pneumoniae (Kumarasamy et al., 2010). Since then, NDM-1 has been detected in numerous countries across all continents, carried by various bacterial species, including Escherichia coli, Klebsiella pneumoniae, Acinetobacter baumanii, and Enterobacter cloacae (Linciano et al., 2019). The gene encoding NDM-1 is often located on plasmids, which are mobile genetic elements that facilitate horizontal gene transfer between bacteria, further accelerating the spread of resistance. This plasmid-mediated transfer allows NDM-1 to disseminate rapidly across different environments, from hospitals to community settings, complicating infection control efforts. The high prevalence of NDM-1 in both clinical and environmental isolates indicates its widespread dissemination and the urgent need for coordinated surveillance and containment strategies.

Due to the serious consequences of NDM-1-mediated antibiotic resistance, it is urgently necessary to find new inhibitors that can effectively target this enzyme. Current treatment options are limited and often involve the use of antibiotics that are more toxic or less effective, such as colistin and tigecycline. These alternatives may have adverse side effects and are not always successful in treating infections caused by NDM-1-producing bacteria. Research efforts are focused on developing new molecules that can inhibit the activity of NDM-1, thereby restoring the efficacy of beta-lactam antibiotics. Strategies include designing small molecules that can chelate the zinc ions essential for NDM-1 activity or developing novel compounds that can bind to the enzyme's active site. While NDM-1 has been proposed to be repressed by chemicals such as captoprils, sulfonamides, boric acid derivatives and thiol-containing small molecules, none of these potential inhibitors have received approval for clinical use. Hence, in order to counteract bacterial antibiotic resistance, it is essential to design/identify, evaluate and further develop novel and potent NDM-1 inhibitors (Klingler et al., 2015).

Computer-assisted drug design (CADD) has revolutionized the field of pharmaceutical research, offering significant advantages in the drug discovery process. However, it also comes with several limitations that need to be addressed for optimal application. These include the complexity and dynamic nature of biological systems, which are often oversimplified in computational models, leading to inaccuracies. The accuracy of predictive models depends on the quality of data, which can be incomplete or biased, and the need for substantial computational resources can be a barrier. Modeling the flexibility of biological targets is challenging, and integrating computational predictions with experimental data is time-consuming and costly. Scalability issues, outdated software and algorithms, and difficulties in accurately predicting ADMET properties further complicate the drug discovery process. Continuous advancements and improved integration with experimental data are needed to fully realize the potential of CADD.

Given the substantial role of NDM-1 in the dissemination of antibiotic resistance, this in silico study is dedicated to identifying potential small molecule inhibitors of NDM-1.

2 Methodology

3 Results

3.1

3.1 Protein structures alignment

After retrieval of protein crystal structures of NDM-1 from PDB, they were aligned on each other to analyze the binding pocket residues. The RMSD difference was less than 0.5 which showed that there were not many deviations in the structures (Fig. 1A). The orientation of the residues of binding pocket was deeply analyzed (Fig. 1B). It was observed that the conformations of Val73, Trp93, His122, Gln123, Asp124, and Lys211 had changed in 6TTC structure.

Close

Fig. 1

The superimposed residues of the NDM-1 crystal structures used for pharmacophore model generation. (A) The alignment of the co-crystal structures, (B) The aligned residues of the binding pocket.

Export to PPT

3.2

3.2 Pharmacophore creation

Binding pocket residues of NDM-1 protein were chosen in order to create a seven-featured pharmacophore model of the receptor cavity. The coordinates in the protein structure for the pharmacophore hypothesis are depicted in Fig. 2A and summarized in Sup. Table S1. Fig. 2B depicts the features of the pharmacophore hypothesis, including the binding pocket cavity.

Close

Fig. 2

Pharmacophore model of NDM-1 for the binding pocket (A) and the cavity of the pocket (B).

Export to PPT

Table 1 Binding free energies for the two ligands complexed to NDM-1 as evaluated by MM-GBSA analyses.

Energy components (kcal/mol)	STK115225-complex	STK107343-complex
ΔEvdW	−52.86 ± 0.34	−59.85 ± 0.37
ΔEele	−2.96 ± 0.25	−8.90 ± 0.36
ΔEGB	13.40 ± 0.22	20.85 ± 0.32
ΔEsurf	−5.50 ± 0.03	−5.63 ± 0.04
ΔGgas	−55.83 ± 0.37	−68.75 ± 0.54
ΔGsolv	7.89 ± 0.21	15.21 ± 0.32
ΔGtotal	−47.93 ± 0.31	−53.54 ± 0.43

3.3

3.3 Virtual screening

Virtual screening of the selected ligands was conducted employing the pharmacophore model of NDM-1. A compound was designated as a hit only when a minimum of four features were found to be matched. Parameters such as vector alignments, RMSD site matching and volume scores were used to rank the final hits from screening in order of their phase fitness scores. Selected ligands’ vector scores ranged from −1.0 to + 1.0., and higher values in the vector score corresponded to greater alignment. On the other hand, volume scores ranged between 0.0 and 1.0, with higher volume scores indicating higher commonality of the volumes of the aligned and reference ligands. Phase screen score threshold was set at 1.4 in order to identify possible hits. The ligands which were selected along with their phase screen score and alignment on the pharmacophore hypothesis are given in Sup. Table S2.

Table 2 Glide g-score scores of the selected ligands.

S. No.	Compound ID	Structure of ligand	Glide score
1.	STK115225		−7.291
2.	STK107343		−7.280
3.	STK215880		−7.262
4.	STK115575		−7.043

3.4

3.4 Molecular docking

The ligands were docked to NDM-1 protein in accordance with the standard methodology of the glide tool. Analyses of molecular interactions between the compounds docked to the receptor were evaluated in consideration with the glide g-score obtained. A cutoff binding energy score of −7 kcal/mol was fixed in order to select the hits. This resulted in the selection of 4 hits with appropriate glide g-score and binding energies (Table 2). Upon analyses of the molecular interactions between the protein and docked hits, we observed that STK115225 made one π-π interaction with His122, and four π-alkyl interactions with Leu65, Val73, Trp93, and His250. STK107343 was found to make one π-sulfur bond with Trp93, one π-σ bond with His122, a π-π bond with His250, and an alkyl bond with Ala215. STK215880 was involved in one π-π interaction with His250, an alkyl bond with Ala215 and a metal acceptor with Zn as well. Lastly, STK115575 elicited three hydrophobic interactions with Ala74, Trp93, and His122. Fig. 3 depicts the ligand-receptor molecular interactions for all the four hits against NDM-1. Analyses of binding modes for the selected ligands docked to NDM-1 were also conducted and are represented in Fig. 4.

Close

Fig. 3

Molecular interactions underlying binding of the four selected ligands (STK115225, STK107343, STK215880 and STK115575) to NDM-1. Magenta, purple and orange lines depict hydrophobic interactions, pi-Sigma bonds, and pi-sulfur interactions, respectively. (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)

Export to PPT

Close

Fig. 4

Probable binding modes of the four selected ligands (STK115225, STK107343, STK215880 and STK115575) represented as sticks in the binding pocket of NDM-1.

Export to PPT

3.5

3.5 MDS analysis

The docked poses of the ligands were superposed on the co-crystal ligand and two compounds were chosen for protein–ligand stability analyses. The aligned poses of the selected compounds are shown in Fig. 5 A-B.

Close

Fig. 5

Probable binding modes for ligands, STK115225 (A; magenta sticks) and STK107343 (B; blue sticks) aligned with NDM-1 co-crystal structure. (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)

Export to PPT

3.5.1

3.5.1 Root mean square deviation (RMSD)

To evaluate the stability of the ligand–protein complexes, RMSD of the alpha carbon atoms in the proteins and ligands was computed from the trajectories (Sargsyan et al., 2017). As depicted in Fig. 6A, complexes for STK115225 and STK107343 and were equilibrated for 10 ns. RMSD of STK115225 complex showed negligible deviations till 20 ns and achieved stability at ∼ 1.25–1.50 Å until the culmination of the MDS. There were minor deviations at the frames generated at 50 and 60 ns. For the STK107343 complex, RMSD attained a value of ∼ 1.75 Å and remained in this range until the culmination of the MDS at 100 ns.

Close

Fig. 6

(A) RMSD plots of ligand–protein complex for the duration of MDS (100 ns). (B) RSMF value plots depicting the residue fluctuations during the simulation period. Cyan and magenta plots represent complexes of STK107343 and STK115225, respectively for both (A) and (B). (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)

Export to PPT

3.5.1.1

3.5.1.1 Root mean square fluctuations (RMSF)

RMSF values are a good measure of the dynamic behavior of amino acid restudies in receptor proteins, when the latter are complexed with their ligands (Martínez, 2015). With the exception of the loop section residues which extended to the maximum value of ∼ 2.5 Å, RMSF of other residues were found to vary by less than 1 Å during the course of the simulation (Fig. 6B). Overall, the results indicate that the protein residues are more-or-less rigid and do not show any major fluctuations. This suggests that the ligands form a stable complex with the protein receptor DNM1.

3.5.1.2

3.5.1.2 Protein-Ligand contacts

Based on MDS analysis, interactions between the ligands and receptor protein, NDM-1 mainly relied on hydrophobic, ionic, and hydrogen bond interactions. Amino acid residues involved in ionic bonds in STK115225-NDM-1 complex were Asp124, Cys208, Lys211, and His250 (Fig. 7A). On the other hand, His120, His122, Asp124, His189, Cys208, Lys211 and His250 were found to make ionic bonds with STK107343 compound (Fig. 7B).

Close

Fig. 7

Ligand-protein interaction during MDS for STK115225 (A) and STK107343 (B). The interacting residues in NDM-1 are displayed as stacked bars.

Export to PPT

3.5.1.3

3.5.1.3 Principal component analysis

PCA was employed to evaluate the dynamic behavior of the protein in both complexes (Fig. 8 A-B). This method aids in determining the collective motion of dynamics trajectories. Plotting the variance proportion against the Eigen values, which depict dynamic motions, was performed. Three principal components (PC) which cover the major fluctuations were plotted. The PCA analysis can demonstrate conformational changes in all groups through basic grouping in PC subspace (Khan et al., 2021). For the STK115225 complex, overall flexibility was recorded at 31.55 %. PC-1 elicited the highest variation of 15.35 %, while PC-2 and 3 had variabilities of 9.04 % and 7.16 %, respectively (Fig. 8A). With regards to the STK107343-complex, PC-1 also showed the maximum variation at 21.68 %. Corresponding variations of 6.8 %, and 5.71 % were observed for PC-2 and 3, respectively. Overall flexibility for the STK107343 complex was found to be 34.19 % (Fig. 8B).

Close

Fig. 8

PCA of ligands STK115225 (A) and STK107343 (B) bound to NDM-1. Blue, white and red represents greatest, intermediate and least amount of flexibility movements. (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)

Export to PPT

3.5.1.4

3.5.1.4 Cross correlation

Correlations amongst NDM-1 protein residues were assessed using a dynamic cross correlation matrix (DCCM; Fig. 9A-B). Most residues elicited positive correlations, but few and short-lived anti-correlated residues were also present. Overall, cross correlation matrix suggests that NDM-1 residues are highly correlated to in presence of the ligands during the MDS.

Close

Fig. 9

Dynamic cross correlation matrix (DCCM) of STK115225 (A) and STK107343 (B) complexed to NDM-1 protein. Cyan represents positively correlated residues, and magenta depicts the anti-correlated ones. Cyan diagonal lines in DCCM indicate positive correlations between topologically proximate residues. (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)

Export to PPT

3.5.1.5

3.5.1.5 MM-GBSA calculations

MM-GBSA analyses were performed in order to calculate the total binding free energy (ΔG_total) of the two ligands, STK115225 and STK107343 to NDM-1 protein. ΔG_total is a good measure of the stability of protein–ligand complex (Du et al., 2011), and lower values indicate enhanced stability of the complex. ΔG_total calculated by MM-GBSA model is dependent on the contributions of the various kinds of ligand–protein interaction energies, including van der Waals (ΔE_vdW) and electrostatic (ΔE_ele) energies, and electrostatic contribution by generalized born (ΔG_GB). ΔG_total for the two complexes are depicted in Table 1.

The ΔE_vdW contribution of STK107343 complex was more than STK115225 complex, and this was same for the electrostatic contribution. Similarly, generalized born contribution was higher for STK107343 when compared to the STK115225 complex. ΔG_total was computed as −47.93 ± 0.31 kcal/mol for STK115225 complex, and −53.54 ± 0.43 kcal/mol for STK107343 complex. The contribution of each energy component in the calculation of ΔG_total is depicted in Supplementary Fig. S1.

4 Discussion

Antibiotic resistance in pathogenic microbial has emerged as a significant global health issue. NDM-1 and its variants have recently attracted a lot of consideration because of their ability to hydrolyze practically all β-lactam antibiotics. These antibiotics include carbapenems, which have been traditionally thought to be last resort antibiotics. Horizontal gene transfer and enhanced spread of NDM-1 in newer bacterial species has aggravated the problems. In this regard, it is critical to identify and clinically evaluate NDM-1 inhibitors which may retard antibiotic resistance in pathogenic microbes and render current antibiotics usable and beneficial (Ahmed et al., 2021; Faheem et al., 2013; Rahman and Khan, 2020). In recent years, the utilization of investigational, current, and experimental drugs in drug development approaches has become increasingly important. Among them, high-throughput computational screening methods such as virtual screening of massive databases concentrating on different chemical principles, network pharmacology, and molecular docking have taken the lead (Löwer and Proschak, 2011; Thangavel and Albratty, 2022). Virtual screening of natural and synthetic molecules has emerged as a vital strategy for modern drug development (Walters et al., 1998). Virtual screening can be used for computational evaluation of novel small molecule chemicals from existent libraries for their binding effectiveness against target proteins. Chemical libraries are filtered using the docking approach, and compounds are ranked according to how well they bind (Lyne, 2002). In the present study, virtual screening in complementation with molecular docking has been employed in order to identify potent high-affinity binding ligands against NDM-1 protein.

The crystal structures of the NDM-1 protein were retrieved from PDB, and were aligned onto each other with reference to the binding pocket residues. The structures were aligned on each other to analyze the binding pocket residues. The RMSD difference was less than 0.5 which showed that there were not many deviations in the structures. Using the pharmacophore hypothesis, virtual screening of the small molecule ligand library obtained from the Vitas-M Laboratory database was carried out. A threshold phase screen score of 1.4 was utilized to identify possible hits, and 44 different compounds scoring higher than this threshold were chosen. These hits were analyzed for their binding to NDM-1 protein by molecular docking studies. This allowed us to predict the ligands’ affinity, orientation, and interaction when complexed with NDM-1. Docking results were evaluated as determined by glide gscores and molecular interaction parameters of the ligands (Silakari and Singh, 2021). In order to select high affinity binding hits and after setting a cutoff value of −7 kcal/mol, 4 hits (STK115225, STK107343, STK215880, and STK115575) were identified and further evaluated.

The docked poses of the ligands were superposed on the co-crystal and two compounds (STK115225, and STK107343) were selected for the protein–ligand stability analysis. STK115225 has been found to interact with the active site residues of NDM-1 such as His122, and His250. In addition, STK115225 was also involved in the interaction with second-shell residues of NDM-1 like Leu65, Val73, and Trp93. Similarly, STK107343 interacted with the catalytic residue Asp124 of NDM-1, along with other active site residues (His122, and His250), and second-shell residue (Trp93). Protein-ligand interactions usually involve receptor flexibility, which enables selectivity in ligand recognition. It follows that it is imperative to comprehend the conformational states of the ligand–protein complexes with the high affinity and specificity binding sites using techniques such as molecular dynamics simulations (MDS). MDS allows for the building of an ensemble of structures that may be used to compute thermodynamic potentials such as binding free energy with great precision (Salmas et al., 2015; Zimmermann et al., 2013). According to RMSF analyses, NDM-1 protein residues exhibited rigidity and did not exhibit substantial variations throughout the simulation duration, indicating high stability of the ligand–protein complex. Further, PCA revealed a large static number of hydrogen bonds. Additionally, the cross-correlation matrix was used to evaluate the correlation between the protein residues. Our results suggest that when the protein residues were bonded to their corresponding molecules during simulation, they showed a strong correlation with one another. Lastly, total binding free energies of the complexed ligands were calculated, and confirmed high stability of the complexes.

5 Conclusions

NDM-1 deactivates β-lactam antibiotics and leads to resistance against several antibiotics generally used in clinical settings. Hence, the present computational analyses explore promising small molecules with high binding affinities that could be exploited for designing NDM-1 inhibitors. Herein, Vitas-M laboratory database was screened for small molecules capable of binding to NDM-1 by generating a structure-based pharmacophore hypothesis, followed by molecular docking and dynamics simulations. The ligands, STK115225 and STK107343 were further analyzed for the stability and affinity of protein–ligand interactions. While we identified key compounds capable of targeting NDM-1, further experimental (in-vitro/-vivo) validation is required to confirm and extend our current findings (Saeidnia et al., 2016). However, in vitro and in vivo experiments are required to determine if these ligands (STK115225 and STK107343) can inhibit the biological function of the NDM1 protein and thus target antibiotic resistance in bacteria.

Ethical approval: This study is exempted from ethical approval as it is a type of in silico study and there is no involvement of patient/healthy subjects during the entire course of the study

Institutional Review Board Statement: Not applicable.

Informed Consent Statement: Not applicable.

CRediT authorship contribution statement

Shafiul Haque: Writing – review & editing, Writing – original draft, Visualization, Validation, Supervision, Software, Resources, Project administration, Methodology, Investigation, Funding acquisition, Formal analysis, Data curation, Conceptualization. Darin Mansor Mathkor: Writing – review & editing, Writing – original draft, Visualization, Validation, Methodology, Investigation, Funding acquisition, Formal analysis, Data curation, Conceptualization. Ayman K. Johargy: Writing – original draft, Visualization, Validation, Methodology, Investigation, Formal analysis, Data curation, Conceptualization. Hani Faidah: Writing – original draft, Visualization, Validation, Methodology, Investigation, Formal analysis, Data curation, Conceptualization. Ahmad O. Babalghith: Writing – original draft, Visualization, Validation, Methodology, Investigation, Formal analysis, Data curation, Conceptualization. Sumyya H. Hariri: Writing – original draft, Visualization, Validation, Methodology, Investigation, Formal analysis, Data curation, Conceptualization. Naif A. Jalal: Writing – original draft, Visualization, Validation, Methodology, Investigation, Formal analysis, Data curation, Conceptualization. Faraz Ahmad: Writing – review & editing, Writing – original draft, Visualization, Validation, Methodology, Investigation, Formal analysis, Data curation, Conceptualization. Farkad Bantun: Writing – original draft, Visualization, Validation, Methodology, Investigation, Formal analysis, Data curation, Conceptualization.