Research Article
Print
Research Article
Computational quest to inhibit enoyl-acyl carrier protein reductase (InhA) enzyme in Mycobacterium tuberculosis from Lannea coromandelica (Houtt.) Merr. metabolite via molecular docking and dynamics
expand article infoMuhammad Askar‡
‡ Politeknik Kesehatan Kementerian Kesehatan Makassar, Makassar, Indonesia

Corresponding author: Muhammad Askar ( muhammadaskar2828@gmail.com )

Academic editor: Emilio Mateev
© 2024 Muhammad Askar.
This is an open access article distributed under the terms of the Creative Commons Attribution License (CC BY 4.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original author and source are credited.
Citation: Askar M (2024) Computational quest to inhibit enoyl-acyl carrier protein reductase (InhA) enzyme in Mycobacterium tuberculosis from Lannea coromandelica (Houtt.) Merr. metabolite via molecular docking and dynamics. Pharmacia 71: 1-10. https://doi.org/10.3897/pharmacia.71.e129151
Open Access

Abstract

This research is a computational exploration to look for natural compounds that can inhibit enoyl-acyl carrier protein reductase (InhA) enzyme from Mycobacterium tuberculosis (MTB). One of Indonesia’s native plants that has been reported to inhibit InhA from MTB is the Lannea coromandelica (Houtt.) Merr., known to contain various active metabolites. However, the molecular activity of the metabolites has not been determined. The aim of this research is the discovery and testing computationally of the binding of metabolites from Lannea coromandelica (Houtt.) Merr. The metabolite ligands are obtained from the natural compound database KNApSAcK, and the 3D structure of the receptor is obtained from the PDB website with the code 1BVR. Subsequently, molecular docking is performed using MGL Tools v.1.5.6 and AutoDock Tools v.4.2.3 software. High-performance computers are used for molecular dynamics simulations with the Gromacs 2016.3 software for a duration of 100 nanoseconds (ns). The docking simulation results show that metabolites show a negative binding energy and close to the value of the native ligand. Moreover, molecular dynamics simulation analysis revealed significant stability of the C1436 ((2R,3R)-3-Hydroxy-5,7,4’-trimethoxyflavanone)-InhA complex over 100 nanoseconds. Molecular dynamics simulations demonstrate that from MM-PBSA value compound C1436 showed MMPBSA binding energy −104.995 kJ/mol, closely approaching the value of the native ligand, which is −142.999 kJ/mol. Furthermore, from the molecular dynamics simulation analysis, C1436 compound demonstrates stability very similar to the native ligand, as observed from the RMSD, RMSF, Rg, RDF, SASA, and PCA analysis. In conclusion, the compound from Lannea coromandelica (Houtt.) Merr. has the potential to serve as a lead compound for inhibiting InhA in MTB.

Keywords

Lannea coromandelica (Houtt.) Merr., Enoyl-acyl carrier protein reductase (InhA), Mycobacterium tuberculosis (MTB), KNApSAcK database, Computational approach

Introduction

Mycobacterium tuberculosis (MTB) is responsible for a global infectious disease that poses a significant threat to adult populations, resulting in fatalities (WHO 2020). Roughly one-third of the global population is believed to harbor MTB infections, with approximately 10% of those individuals expected to experience active tuberculosis (TB) at some point in their lives. The conventional approach to TB treatment involves the use of drug combinations to mitigate the risk of drug resistance development. However, inadequate and prolonged treatment plans have contributed to the rise of multidrug-resistant tuberculosis (MDR-TB) and extensively drug-resistant TB (XDR-TB). Addressing MDR-TB and XDR-TB requires prolonged treatment periods, resulting in significant expenses and frequently causing serious side effects. This situation complicates efforts to control tuberculosis. Consequently, there is an imperative requirement for the discovery of novel and potent anti-TB medications capable of eradicating drug-resistant mycobacteria (Singh et al. 2012).

Isoniazid (INH), Ethionamide, and Triclosan, recognized as key anti-TB drugs, have been pinpointed as substances that specifically focus on the enoyl-acyl carrier protein reductase (InhA) enzyme encoded by InhA (Agarwal et al. 2019). InhA plays a crucial role in the synthesis of long-chain fatty acid derivatives, which serve as vital precursors for mycolic acid and are essential for cell wall construction and various physiological functions in Mycobacterium tuberculosis (MTB). This enzyme plays a crucial role in the NADH-dependent reduction of double bonds in 2-trans-enoyl acyl carrier proteins, a pivotal step in the type II fatty acid synthase (FAS-II) system. INH is classified as a prodrug, requiring activation by the MTB catalase-peroxidase enzyme (Kumar and Navatheesh 2017). Nicotinic acid, the activated form of INH, operates by establishing covalent bonds with NADH within the active site of the protein. This process ultimately results in the production of nicotinamide adenine dinucleotide (NAD). It is noteworthy that the FAS-II system in MTB differs significantly from the mammalian fatty acid synthesis pathway, making the InhA enzyme an appealing target for drug discovery endeavors (Lefebvre et al. 2020).

One of the plants found in tropical areas of Indonesia is Javanese wood, with the Latin name Lannea coromandelica (Houtt.) Merr. This plant has antihypertensive, wound healing and antimicrobial properties (Zhou et al. 2021). Ethnopharmacologically by Indian people, Lannea coromandelica (Houtt.) Merr. has been used in the treatment of infectious diseases such as worms and skin infections (Kaur et al. 2013). The antioxidant, antimicrobial and antithrombolytic properties of ethanol extracts of Lannea coromandelica (Houtt.) Merr. leaves (Alam et al. 2017) have been reported. There is potential activity of Lannea coromandelica (Houtt.) Merr. which can inhibit the growth or kill microorganisms. Previous research results reported that, the Lannea coromandelica (Houtt.) Merr. clicka has inhibitory activity against the growth of MTB. The extracts and partitions obtained from these plants have been tested for anti-TB activity using the method of detection and separation (MODS) and show that methanol extract, hexane soluble and hexane insoluble partitions can inhibit the growth of MTB (Hamzah et al. 2017; Widawati 2022). However, until now there has been no research that shows the mechanism of action of the bioactive compounds in Lannea coromandelica (Houtt.) Merr., and their potential for the specific receptor of MTB, namely InhA, through a computational chemistry approach.

The objective of this research is to discover novel and promising lead compounds that can serve as inhibitors of InhA, a vital component of the FAS-II MTB system (Soeroto et al. 2019). To accomplish this objective, a comprehensive strategy incorporating high-throughput virtual docking screening (VS) and dynamic simulation studies was utilized to explore analogs of natural compounds. This comprehensive approach aimed to identify potential InhA inhibitors and advance their development in the field. Through the utilization of computational techniques such as molecular docking and molecular dynamics simulations, this research aimed to assess the binding affinity and stability of analog natural compounds with InhA. This approach aimed to provide insights into the potential of these compounds as inhibitors of this crucial enzyme (Tierno et al. 2023). The discovery of new InhA inhibitors holds great promise in combating tuberculosis, especially drug-resistant strains, by offering alternative therapeutic avenues and addressing the growing global health concern posed by this infectious disease. This study adds to the continual endeavors in discovering novel solutions for tuberculosis treatment and underscores the importance of computational approaches in the process of drug discovery and development.

Materials and methods

Preparation of macromolecule protein structure

In this study, we opted for the crystallographic structure of enoyl-acyl carrier protein reductase (InhA) enzyme from the Protein Data Bank (PDB ID: 1BVR), based on our prior research findings (Rozwarski et al. 1999). The main objective of this inquiry was to establish a detailed structure of the target protein, with a specific focus on the inhibitor in its bound state, particularly highlighting nicotinamide adenine dinucleotide (NAD) as the native ligand. To ensure the utmost quality and reliability, the crystallographic structure underwent a thorough cleaning process, involving the removal of bound ligands, water molecules, and chain B. Following this, we meticulously added polar hydrogen atoms and Kollman charges, employing techniques available in MGL Tools v.1.5.6 and AutoDock Tools v.4.2.3 (Morris et al. 2010; Huey et al. 2012).

Preparation of natural compound library

The chosen ligands were derived from the terrestrial natural products of Lannea coromandelica (Houtt.) Merr. The structures of these compounds were acquired from the KnpaSack database (http://www.knapsackfamily.com/KNApSAcK/ (Afendi et al. 2012) and http://knapsack3d.sakura.ne.jp/) (Nakamura et al. 2013).

Molecular docking-based virtual screening

The Mycobacterium tuberculosis (MTB) receptor chosen was the crystallographic structure of enoyl-acyl carrier protein reductase (InhA) enzyme from the Protein Data Bank (PDB ID: 1BVR) based on high resolution and bound inhibitor (nicotinamide adenine dinucleotide (NAD)) as native ligand. Bound ligands, water molecules, and B chains are removed (Biela et al. 2012). The approach involved introducing polar hydrogen atoms and Kollman charges using the procedure implemented in MGL Tools v.1.5.6 and AutoDock Tools v.4.2.3. The protein target’s active site was defined by the coordinates: center X = 19.579, center Y = 9.067, center Z = 14.345, corresponding to the binding location of the native ligand. The box dimensions were set to 64 × 60 × 60 Å, with a grid spacing of 0.375 Å. After aligning the X-ray conformation of NAD with the best conformation acquired through the redocking process, the RMSD value for InhA’s grid box parameters was 1.8 Å (Nurisyah et al. 2024).

Identification of the protein-ligand complex

The visualization and analysis of protein-ligand complexes were carried out using software tools such as Pymol (version 1.7.4.5 Edu), the Molecular Graphics System, and Discovery Studio Visualizer (D.S.V.) version 16.1.0.153500 (Baroroh et al. 2023). In these software applications, an examination of interactions, encompassing hydrogen bonds, polar interactions, and hydrophobic interactions between the proteins and ligands, was observed and illustrated. These software platforms enabled a comprehensive assessment of the binding modes and conformational changes of the ligands within the protein’s active site, shedding light on their potential as inhibitors (Rizkita et al. 2024). This structural analysis is pivotal in understanding the mechanisms of ligand-protein interactions and aids in the rational design of novel inhibitors with enhanced binding affinity and specificity.

Molecular Dynamics (MD) simulations

The Gromacs 2016.3 software was employed for the molecular dynamic simulation of the protein-ligand complex, and this computational study utilized various tools and techniques to ensure accuracy (Abraham et al. 2015). Initial topology and ligand parameters were generated using the Antechamber Python Parser Interface (ACPYPE) tools (Sousa da Silva and Vranken 2012), while both the ligand and protein were assigned AMBER99SB-ILDN force fields (Vermaas et al. 2016). To account for electrostatic forces, the Ewald Particle Mesh method was employed, ensuring precise calculations (Isele-Holder et al. 2012). The system was neutralized by the addition of Na+ and Cl- ions, and the solvation process utilized the TIP3P water cube model (Smith et al. 2015). The simulation procedure involved various essential stages, beginning with an initial minimization phase that included 10,000 steps. This phase utilized the steepest descent method for the first 5,000 steps, followed by the conjugate gradients method for the subsequent 5,000 steps (Parikesit and Tambunan 2018). Subsequently, the entire system was subjected to a heating phase, raising the temperature to 310 K over 100 ps, starting from an initial temperature of 0 K (Baral et al. 2021). This temperature regulation was achieved using Langevin dynamics (Hicks et al. 2021). The primary production phase of the simulation took place under constant temperature and pressure conditions, with a temperature held at 310 K and pressure at 1 atm. This phase employed a time step of 2 fs, and the SHAKE algorithm was applied to constrain all hydrogen bonds. Throughout the 100 ns (100,000 ps) simulation, trajectory coordinates were saved at intervals of 20 ps. To comprehensively assess the molecular dynamics simulation results, various analyses were performed, including the determination of root mean squared deviation (RMSD), root mean squared fluctuation (RMSF), radius of gyration (Rg), radial distribution function (RDF), surface area solvent accessible (SASA), principal component analysis (PCA), examination of secondary structure changes, and hydrogen bond analysis. These analyses provide valuable insights into the stability and dynamic behavior of the protein-ligand complex during the simulation (Fakih et al. 2024).

Binding free energy (MM-PBSA) calculation

The Molecular Mechanics-Poisson-Boltzmann Surface Area (MM/PBSA) calculations were carried out using the g_mmpbsa package, which is integrated into the Gromacs 2016.3 software (Kumari et al. 2014). Energy calculations were performed by averaging the data from 100 snapshots, each extracted at ten ps intervals, encompassing the MD trajectory spanning from 0 to 100 ns. To compute the polar desolvation energy, the Poisson–Boltzmann equation was employed with a grid size set at 0.5 Å (Genheden and Ryde 2015). In this calculation, the dielectric constant for the solvent was configured to 80, representing water as the solvent. The nonpolar contribution to the energy was determined by assessing the solvent-accessible surface area with a solvent radius of 1.4 Å. The comprehensive analysis of protein-ligand free energy was carried out based on the results obtained from the molecular dynamics simulations, specifically when the interactions within the complex achieved a stable state (Ren et al. 2020).

Results

Molecular docking-based virtual screening

In this study, various compounds have been virtually screened to assess their ability to bind to the InhA receptor, which has a key role in the tuberculosis infection process. This step is also carried out to obtain the initial conformation used in molecular dynamics simulations. The docking results showed that morin showed better (lower) binding energy compared to native ligand, as well as several other metabolites contained in Lannea coromandelica (Houtt.) Merr. Where morin shows a bond energy of −9.79 kcal/mol while native is −9.33 kcal/mol (Table 1). However, other metabolites also show negative binding energies and are close to the values of the native ligand. It can be interpreted that the metabolites has the potential to be an effective inhibitory agent in interfering with InhA activity on tuberculosis bacteria. Therefore, the conformations based on docking of the three metabolites were chosen to be analyzed using molecular dynamics and compared with the native ligand.

Table 1.
Download as
CSV
XLSX

Binding energy from molecular docking-based virtual screening results of metabolite in Lannea coromandelica (Houtt.) Merr.

Compound name Knapsack ID Binding energy
Native − −9.33 kcal/mol
(2R,3R)-3-Hydroxy-5,7,4’-trimethoxyflavanone C1436 −7.90 kcal/mol
Morin C4624 −9.79 kcal/mol
(2R,3S)-3,5,3’-Trihydroxy-7,4’-dimethoxyflavanone C14368 −7.67 kcal/mol

Structural exploration from Molecular Dynamics (MD) systems

Root-Mean-Square Deviation (RMSD) and Root-Mean-Square Fluctuation (RMSF)

The RMSD plots were analyzed to identify changes in the conformation and structure of the InhA complex and its ligands (NAD, C1436, C4624, and C14368) throughout the 100 ns MD simulations. During the initial 0–50 ns simulation period, all complexes demonstrated similar fluctuations. However, beyond 50–100 ns, the InhA-C1436 complex exhibited more pronounced changes compared to the InhA-NAD, InhA-C4624, and InhA-C14368 complexes. Conversely, the InhA-C4624 and InhA-C14368 complexes showed instability but remained comparable to the InhA-NAD complex. Notably, the InhA-C4624 and InhA-C14368 complexes demonstrated lower fluctuations than the InhA-NAD complex towards the end of the simulation (50–100 ns). The mean RMSD values for the InhA-NAD, InhA-C1436, InhA-C4624, and InhA-C14368 complexes were 0.275 nm, 0.321 nm, 0.262 nm, and 0.273 nm, respectively. These results suggest that the average RMSD value between the identified ligands and NAD does not exhibit significant differences (Fig. 1).

Figure 1. 

The examination of RMSD and RMSF plots during the MD simulations includes representations of complexes with the protein’s native ligand (NAD) and the protein’s hit ligands, each distinguished by different colors. NAD is represented in black, C1436 in red, C4624 in yellow, and C14368 in green.

RMSF plots offer insights into the flexibility of the protein over the trajectory. These plots revealed that all complexes displayed comparable fluctuations, especially in the region of amino acids situated within the active site of the InhA protein. This includes Gly14, Ser20, Leu63, Val65, Ser94, Gly96, Lys165, Ile194, and Thr196, which establish hydrogen bonds with the ligands. Similar fluctuation patterns were observed in amino acid residues within the sequence 190–210, also situated in the active site of InhA. Nevertheless, variations in fluctuations were noticeable in amino acid residues that do not directly interact with the ligand. This was observed in amino acid sequences 48 (InhA-C1436 and InhA-CC14368 complexes), 112 (InhA-C1436 and InhA-CC14368 complexes), 124 (InhA-C1436 complexes, InhA-C4624, and InhA-CC14368), 152 (InhA-C4624 complex), and 176 (InhA-C1436 and InhA-CC14368 complex).

Radius gyration (Rg) and Solvent Accessible Surface Area (SASA)

The purpose of assessing the radius of gyration (Rg) is to evaluate the stability of the complex, specifically whether it remains stable in a folded or open conformation throughout the simulation. The Rg plot exhibits analogous characteristics among the protein-native ligand complex (InhA-NAD) and the protein-bound ligands (InhA-C1436, InhA-C4624, and InhA-C14368). The average Rg values for each ligand are 1.818 nm, 1.809 nm, 1.819 nm, and 1.817 nm. The Rg values for the protein-bound ligand complexes display remarkable resemblances and exhibit stable folded structures when compared to the protein-native ligand complexes (Fig. 2). Rg values below 2 nm indicate that the molecule or protein-ligand complex tends to maintain a more compact and compressed conformation or structure during molecular dynamics simulations. This suggests that the interaction between the protein and the ligand likely results in a stable complex with minimal conformational changes throughout the observed simulation period.

Figure 2. 

The analysis of Rg and SASA plots during the course of MD simulations is illustrated. Different colors are employed to represent the complexes, with the protein’s native ligand (NAD) shown in black, C1436 in red, C4624 in yellow, and C14368 in green.

Moreover, a evaluation of solvent-accessible surface area (SASA) was performed to gauge the extent to which protein conformational changes, occurring during the simulation, are exposed to water molecules. The SASA plots produced similar results for both complexes containing the native protein ligand and those with the hit protein ligands. The average SASA values were 129.126 nm^2, 129.379 nm^2, 127.542 nm^2, and 127.820 nm^2 for the native protein ligand (InhA-NAD) and the hit protein ligands (InhA-C1436, InhA-C4624, and InhA-C14368), respectively (Fig. 2). SASA values exceeding 125 nm^2, as observed in the molecular dynamics simulations, suggest that the surface of this protein-ligand complex tends to engage in more interactions with the solvent and its surrounding environment throughout the simulation. This may imply a higher degree of dynamic conformational changes and increased exposure to the solvent during the simulation.

Principal Component Analysis (PCA)

Principal component analysis (PCA) was conducted to identify noteworthy fluctuations within the protein-ligand complex. The extent and direction of the eigenvectors governing these motions were evaluated by generating a hedgehog plot, offering a visual depiction of motion throughout the trajectory. In this study, all eigenvectors were employed to calculate motion. The eigenvectors indicated percentages of 36.91%, 25.76%, 88.44%, and 24.95% for the protein-native ligand complex (InhA-NAD), the protein-hit ligand complexes (InhA-C1436, InhA-C4624, and InhA-C14368), and the first ligand of the 20 ns simulation, respectively. Moreover, the protein-ligand complex exposed to the ligand formed a stable complex with reduced motion, closely resembling the motion observed in the reference InhA-NAD complex. The investigation also encompassed an assessment of complex dynamics through 2D projection trajectory plotting. A smaller cluster footprint indicates a more stable complex, while larger clusters suggest a less stable protein-ligand complex. The protein-exposed ligand complex was found to occupy less space compared to the protein’s native ligand complex. Nonetheless, the patterns exhibited tend to be similar (Fig. 3).

Figure 3. 

Assessment of eigenvalues from the covariance matrix and the 2D projection of trajectory data during MD simulations was performed. The protein-ligand complexes, including the protein-native ligand (NAD) and protein-hit ligands, were depicted using distinct colors: NAD in black, C1436 in red, C4624 in yellow, and C14368 in green.

Radial Distribution Function (RDF)

The Radial Distribution Function (RDF) is employed to depict the relative arrangement of atom pairs within a molecular system. It serves as a valuable metric for analyzing molecular dynamics simulations and various chemical systems. The RDF quantifies how frequently one atom comes into contact with another atom within a specified distance range. In molecular simulations, RDF is frequently used to assess interactions among specific molecules in a system, such as the interactions between water molecules and solute molecules in a solution. As observed in the earlier analysis, there was no notable disparity in RDF values between the protein-native ligand complex (InhA-NAD) and the protein-hit ligand complex (InhA-C1436, InhA-C4624, and InhA-C14368), with respective values of 5,102 g(r), 5,099 g(r), 4,799 g(r), and 4,916 g(r) (Fig. 4). When RDF values from molecular dynamics simulations exhibit no significant differences, it indicates that the spatial distribution of atoms or molecules within the system maintains similar characteristics and undergoes minimal alterations throughout a specified simulation period.

Figure 4. 

The assessment of RDF was carried out during the course of molecular dynamics (MD) simulations. Different colors were used to represent the protein-ligand complexes, which included the native ligand (NAD) and the hit ligands: NAD in black, C1436 in red, C4624 in yellow, and C14368 in green.

Binding Free Energy (MM-PBSA) calculation InhA-ligands

The MM-PBSA energy calculations were conducted to predict the individual energy contributions to the stability of the ligand-protein complex. These calculations were performed based on the data obtained from 100 ns of molecular dynamics simulations. The results indicate that the InhA-C1436 and InhA-C14368 complexes exhibit favorable predicted binding affinities. Specifically, C1436 demonstrated the lowest binding free energy at −104.995 kJ/mol, while C14368 exhibited a binding free energy of −95.528 kJ/mol. However, it is noteworthy that both the InhA-C1436 and InhA-C14368 complexes did not exhibit a superior predicted binding affinity and stability compared to NAD, the original ligand, which possessed a binding free energy of −142.999 kJ/mol. Moreover, the polar solvation energy was less favorable to the overall energy of the complex systems. In contrast, van der Waals energy and electrostatic energy played significant roles in contributing to the negative binding free energy of the complexes during the molecular dynamics simulations (Table 1). On the other hand, the InhA-C4624 complex demonstrated a positive binding free energy of 38.580 kJ/mol. This positive energy arises from the interaction between the C4624 ligand and the InhA protein, which is less favorable in terms of free binding.

Hydrogen bonds formed between the receptor and ligands

The analysis of hydrogen bond occupancy in the InhA ligand complex was conducted by assessing the percentage of time hydrogen bonds were formed during a 100 ns simulation. Hydrogen bonds with an occupancy percentage exceeding 50% were considered for evaluation. Among the ligands, NAD exhibited the highest number of hydrogen bonds, with a total of 28 (Fig. 5). These hydrogen bonds involved various amino acid residues, including Ser20 (22.59%), Ile21 (11.28%), Ala22 (4.57%), Arg43 (1.12%), Val65 (0.20%), Gln66 (0.12%), Ser94 (0.60%), Gly96 (23.36%), Lys165 (10.99%), Ile194 (7.96%), Thr196 (15.93%), Leu197 (0.24%), And Ala198 (1.43%) As Donors, As Well As Gly40 (0.01%), Phe41 (1.78%), Asp42 (0.53%), Leu63 (1.51%), Asp64 (0.69%), Gln66 (4.02%), Ser94 (33.46%), Ile95 (0.01%), Gly96 (15.23%), Met147 (6.80%), Asp148 (0.49%), Gly192 (2.74%), Ile194 (19.13%), Leu197 (0.51%), and Met199 (0.01%) as acceptors.

Figure 5. 

The assessment of hydrogen bonds was performed during the entire course of molecular dynamics (MD) simulations. The protein-ligand complexes, consisting of the native ligand (NAD) and the hit ligands, were visually distinguished using distinct colors: NAD represented in black, C1436 in red, C4624 in yellow, and C14368 in green.

Conversely, during the simulation, C1436, C4624, and C14368 established 8, 18, and 9 hydrogen bonds, respectively, with specific amino acid residues. C1436 displayed robust hydrogen bonds with amino acid residues Ile21 (0.03%), Phe41 (0.10%), Ser94 (0.09%), Gly96 (45.82%), Tyr158 (1.05%), Lys165 (0.01%), and Ala198 (0.32%). On the other hand, C4624 formed substantial hydrogen bonds with amino acid residues Gly14 (2.90%), Thr39 (13.96%), Gly40 (0.10%), Asp42 (0.57%), Arg43 (28.49%), Leu63 (3.62%), Val65 (0.69), Ser94 (0.12), Gly96 (3.29%), Gln100 (0.29%), Lys165 (4.46%), Gly192 (0.09%), and Ile194 (0.16%). Similarly, C14368 exhibited significant hydrogen bonds with amino acid residues Gly14 (2.97%), Ile15 (2.60%), Thr39 (40.19%), Arg43 (0.07%), Ser94 (0.07%), Gly96 (18.57%), and Ala198 (0.08%).

Visualization of conformational changes

The results showed that the native ligand underwent a conformational change when interacting with the InhA Mycobacterium tuberculosis (MTB) TB receptor from the beginning of the simulation at 0 ns (Fig. 6). This change is in the form of a shift or rotation in the purine group of the molecule as seen in the images captured at 0 ns to 50 and 100 ns of simulation (Fig. 6A). It can be seen that in compound C1436, the benzene methyl ether group changes slightly from its initial position, although it appears that the change that occurs is not up to 90° from its initial position, and then continues to be constant until 100 ns of simulation.

Figure 6. 

Visualization of conformational changes of native ligand (A) and C1436 compound (B) during the molecular dynamics simulation observed at 0 ns, 50 ns, and 100 ns.

Table 2.
Download as
CSV
XLSX

Energy summation results of MM-PBSA calculations for the native ligand (NAD) and analogs of natural compounds.

Molecule Van der Waals Energy (kJ/mol) Electrostatic Energy (kJ/mol) Polar Solvation Energy (kJ/mol) S.A.S.A. Energy (kJ/mol) Bond Prediction Energy (kJ/mol)
Native − 314.845 ± 30.822 −201.204 ± 30.014 399.895 ± 44.128 −26.845 ± 1.543 −142.999 ± 25.304
C1436 −176.366 ± 13.287 −33.705 ± 8.934 122.719 ± 13.317 −17.643 ± 0.958 −104.995 ± 12.288
C4624 −308.309 ± 22.068 −41.037 ± 62.171 415.003 ± 75.381 −27.077 ± 1.422 38.580 ± 32.629
C14368 −172.319 ± 14.579 −55.192 ± 14.730 149.042 ± 22.611 −17.060 ± 0.913 −95.528 ± 18.405

Discussion

The utilization of virtual drug screening has proven its efficacy in swiftly identifying potential lead molecules from a compound library. In our study, we conducted virtual drug screening and dynamic simulation investigations aimed at identifying promising lead compounds that hold potential as candidates for tuberculosis (TB) drug development, with the aim of enhancing TB treatment strategies. We have successfully incorporated natural analog compounds into the Mycobacterium tuberculosis enoyl-acyl carrier protein reductase (MTB-InhA). Among these analogs, three compounds exhibited a substantial binding affinity, comparable to the native ligand (NAD), signifying their potential as potent inhibitors. Through computational methodologies, we have demonstrated that the binding of these three natural compound analogs to MTB-InhA induces notable conformational changes in the protein. Our computational docking investigations have revealed that these analog compounds establish strong interactions with crucial amino acid residues, firmly anchoring themselves within the substrate binding site. These vital interacting residues include Tyr158, Ile194, and Thr196. Notably, Thr196, situated within the binding loop between NAD and InhA, plays a pivotal role in stabilizing the active site NAD. Moreover, the hydroxyl group of Ile194 facilitates hydrogen bonding between InhA and NAD. Additionally, Tyr158 assumes a significant role in these interactions, and mutations in Ala94 have been linked to Mycobacterium tuberculosis resistance to Isoniazid (INH), Ethionamide, and Triclosan.

Subsequently, we performed a 100 ns molecular dynamics (MD) simulation on the initially docked complexes to evaluate the binding affinity and flexibility of three compound analogs—specifically, C1436, C4624, and C14368—with MTB-InhA. Our observations indicated that the binding configurations maintained their energetic stability throughout the MD simulation. We assessed the impact of ligand-induced conformational changes on the protein during MD simulation by monitoring fluctuations in amino acid residues. The highest fluctuations were primarily observed in the amino acid chain corresponding to the loop and C-terminal regions of the InhA protein. Remarkably, the three compound analogs exhibited conformational changes resembling those induced by NAD within the InhA protein structure. Notably, C1436 and C14368 exhibited fluctuations highly resembling those of NAD, particularly in the region spanning residues 100 to 115, where other ligands exerted a comparatively minor influence. Furthermore, C1436 and C14368 displayed minimal fluctuations towards the InhA protein, while C4624 notably impacted the C-terminal residues.

The rise of multidrug-resistant Mycobacterium tuberculosis (MTB) has resulted in a significant increase in mortality rates, severely constraining treatment options. Previous research has clarified that MTB-InhA plays a crucial role in the development of resistance to frontline drugs like INH, Ethionamide, and Triclosan, which are essential in tuberculosis (TB) treatment. Notably, around 80% of MTB isolates globally exhibiting phenotypic resistance to Ethionamide and Triclosan showcase mutations either in codon 315 of the katG gene or at position −15 in the inhA promoter. While our prior knowledge attributed drug activation to the catalase-peroxidase enzyme encoded by katG in the case of INH, compound analogs do not require this activation step via the enzyme. Consequently, these compounds may retain their efficacy against MTB resistance stemming from mutations in katG. Furthermore, we hypothesize that the observed reduction in binding affinity of natural compound analogs may result from perturbations in the hydrogen-bonding network within the InhA active sites, particularly affecting Ile194 and conceivably Tyr158. These perturbations could directly impact the compounds’ binding capability, as previously discussed.

Conversely, C1436, C4624, and C14368, the three analogous compounds, established hydrogen bonds with Ser94 and Ile16 of InhA. Another investigation has proposed that mutations S94A and I16T could potentially induce conformational alterations in the enzyme’s active site upon ligand binding, leading to reduced binding affinity. These findings imply that the three most promising analogous compounds may exhibit binding affinity with MTB-InhA, even in the presence of mutations distinct from those previously described. Nevertheless, comprehensive in silico and in vitro studies are imperative to assess the efficacy of these three lead ligands, particularly against various MTB strain mutations that confer resistance to INH, Ethionamide, and Triclosan.

Conclusions

In this study, we underwent computational exploration to search for natural compounds that have the potential to inhibit enoyl-acyl carrier protein reductase (InhA) enzyme from Mycobacterium tuberculosis (MTB), which is an important target in tuberculosis drug development. We used a molecular docking approach and molecular dynamics simulations to examine interactions between natural compounds and InhA. The results of docking simulations show that the C1436 compound ((2R,3R)-3-Hydroxy-5,7,4’-trimethoxyflavanone) stands out with a binding energy close to the native ligand. Molecular dynamics simulation analysis strengthens these findings by demonstrating the stability of the C1436-InhA complex over 100 ns.

Acknowledgments

The author thanks Department of Medical Laboratory Technology, Politeknik Kesehatan Kementerian Kesehatan Makassar for the support that enabled the completion of this research.

This research received no external funding.

References

  • Abraham MJ, Murtola T, Schulz R, Páll S, Smith JC, Hess B, Lindahl E (2015) GROMACS: High performance molecular simulations through multi-level parallelism from laptops to supercomputers. SoftwareX 1–2: 19–25. https://doi.org/10.1016/j.softx.2015.06.001
  • Afendi FM, Okada T, Yamazaki M, Hirai-Morita A, Nakamura Y, Nakamura K, Ikeda S, Takahashi H, Altaf-Ul-Amin M, Darusman LK, Saito K, Kanaya S (2012) KNApSAcK family databases: integrated metabolite-plant species databases for multifaceted plant research. Plant & Cell Physiology 53: e1. https://doi.org/10.1093/pcp/pcr165
  • Agarwal A, Katoch CDS, Kumar M, Dhole TN, Sharma YK (2019) Evaluation of Microscopic observation drug susceptibility (MODS) assay as a rapid, sensitive and inexpensive test for detection of tuberculosis and multidrug resistant tuberculosis. Medical Journal, Armed Forces India 75: 58–64. https://doi.org/10.1016/j.mjafi.2018.03.011
  • Alam MB, Kwon K-R, Lee S-H, Lee S-H (2017) Lannea coromandelica (Houtt.) Merr. Induces Heme Oxygenase 1 (HO-1) Expression and Reduces Oxidative Stress via the p38/c-Jun N-Terminal Kinase–Nuclear Factor Erythroid 2-Related Factor 2 (p38/JNK–NRF2)-Mediated Antioxidant Pathway. International Journal of Molecular Sciences 18(2): 266. https://doi.org/10.3390/ijms18020266
  • Baral K, Adhikari P, Jawad B, Podgornik R, Ching W-Y (2021) Solvent effect on the structure and properties of RGD peptide (1FUV) at body temperature (310 K) using ab initio molecular dynamics. Polymers 13(19): 3434. https://doi.org/10.3390/polym13193434
  • Baroroh U, Muscifa ZS, Destiarani W, Rohmatullah FG, Yusuf M (2023) Molecular interaction analysis and visualization of protein-ligand docking using Biovia Discovery Studio Visualizer. Indonesian Journal of Computational Biology (IJCB) 2: 22–30. https://doi.org/10.24198/ijcb.v2i1.46322
  • Biela A, Khayat M, Tan H, Kong J, Heine A, Hangauer D, Klebe G (2012) Impact of ligand and protein desolvation on ligand binding to the S1 pocket of thrombin. Journal of Molecular Biology 418(5): 350–366. https://doi.org/10.1016/j.jmb.2012.01.054
  • Fakih TM, Ritmaleni R, Zainul R, Muchtaridi M (2024) Molecular docking-based virtual screening and computational investigations of biomolecules (curcumin analogs) as potential lead inhibitors for SARS-CoV-2 papain-like protease. Pharmacia 71: 1–19. https://doi.org/10.3897/pharmacia.71.e123948
  • Hamzah N, Wahid N, Dhuha S, Tahir KA, Febrianti AP, Ismail I (2017) Aktivitas Inhibisi Pertumbuhan Plasmodium falciparum dan Mycobacterium tuberculosis dari Ekstrak dan Partisi Klika Kayu Jawa (Lannea coromandelica [Houtt.] Merr.). Journal of Pharmaceutical and Medicinal Sciences 2: 85–90.
  • Hicks A, MacAinsh M, Zhou H-X (2021) Removing thermostat distortions of protein dynamics in constant-temperature molecular dynamics simulations. Journal of Chemical Theory and Computation 17: 5920–5932. https://doi.org/10.1021/acs.jctc.1c00448
  • Huey R, Morris GM, Forli S (2012) Using AutoDock 4 and AutoDock Vina with AutoDockTools: A Tutorial. The Scripps Research Institute Molecular.
  • Isele-Holder RE, Mitchell W, Ismail AE (2012) Development and application of a particle-particle particle-mesh Ewald method for dispersion interactions. Journal of Chemical Physics 137: 174107. https://doi.org/10.1063/1.4764089
  • Kaur R, Jaiswal ML, Jain V (2013) Protective effect of Lannea coromandelica Houtt. Merrill. against three common pathogens. Journal of Ayurveda and Integrative Medicine 4: 224–228. https://doi.org/10.4103/0975-9476.123706
  • Kumari R, Kumar R, Consortium OSDD, Lynn A (2014) g _ mmpbsa - A GROMACS tool for MM-PBSA and its optimization for high-throughput binding energy calculations. Journal of Chemical Information and Modeling 54(7): 1951–1962. https://doi.org/10.1021/ci500020m
  • Morris GM, Huey R, Lindstrom W, Sanner MF, Belew RK, Goodsell DS, Olson AJ (2010) AutoDock4 and AutoDockTools4: Automated docking with selective receptor flexibility. Reference-36 docking simulation. Journal of Computational Chemistry 30(16): 2785–2791. https://doi.org/10.1002/jcc.21256
  • Nakamura K, Shimura N, Otabe Y, Hirai-Morita A, Nakamura Y, Ono N, Ul-Amin MA, Kanaya S (2013) KNApSAcK-3D: a three-dimensional structure database of plant metabolites. Plant & Cell Physiology 54: e4. https://doi.org/10.1093/pcp/pcs186
  • Nurisyah, Ramadhan DSF, Dewi R, Asikin A, Daswi DR, Adam A, Chaerunnimah, Sunarto, Rafika Artati, Fakih TM (2024) Targeting EGFR allosteric site with marine-natural products of Clathria sp.: A computational approach. Current Research in Structural Biology 7: 100125. https://doi.org/10.1016/j.crstbi.2024.100125
  • Ren J, Yuan X, Li J, Lin S, Yang B, Chen C, Zhao J, Zheng W, Liao H, Yang Z, Qu Z (2020) Assessing the performance of the g_mmpbsa tools to simulate the inhibition of oseltamivir to influenza virus neuraminidase by molecular mechanics Poisson–Boltzmann surface area methods. Journal of the Chinese Chemical Society 67(1): 46–53. https://doi.org/10.1002/jccs.201900148
  • Rizkita AD, Dewi SA, Fakih TM, Lee CC (2024) Effectiveness of sesquiterpene derivatives from Cinnamomum genus in nicotine replacement therapy through blocking acetylcholine nicotinate: a computational analysis. Journal of Biomolecular Structure and Dynamics 1–14. https://doi.org/10.1080/07391102.2024.2305315
  • Rozwarski DA, Vilchèze C, Sugantino M, Bittman R, Sacchettini JC (1999) Crystal structure of the Mycobacterium tuberculosis enoyl-ACP reductase, InhA, in complex with NAD+ and a C16 fatty acyl substrate. Journal of Biological Chemistry 274(22): 15582–15589. https://doi.org/10.1074/jbc.274.22.15582
  • Singh S, Kumar P, Sharma S, Mumbowa F, Martin A, Durier N (2012) Rapid Identification and Drug Susceptibility Testing of Mycobacterium tuberculosis: Standard Operating Procedure for Non-Commercial Assays: Part 1: Microscopic Observation Drug Susceptibility Assay v2.4.12. Journal of Laboratory Physicians 4: 101–111. https://doi.org/10.4103/0974-2727.105592
  • Smith MD, Rao JS, Segelken E, Cruz L (2015) Force-Field Induced Bias in the Structure of Aβ21–30: A Comparison of OPLS, AMBER, CHARMM, and GROMOS Force Fields. Journal of Chemical Information and Modeling 55(12): 2587–2595. https://doi.org/10.1021/acs.jcim.5b00308
  • Soeroto AY, Lestari BW, Santoso P, Chaidir L, Andriyoko B, Alisjahbana B, Van Crevel R, Hill PC (2019) Evaluation of Xpert MTB-RIF guided diagnosis and treatment of rifampicin-resistant tuberculosis in Indonesia: A retrospective cohort study. PLOS ONE 14(2): e0213017. https://doi.org/10.1371/journal.pone.0213017
  • Tierno D, Grassi G, Zanconati F, Bortul M, Scaggiante B (2023) An overview of circulating cell-free nucleic acids in diagnosis and prognosis of triple-negative breast cancer. International Journal of Molecular Sciences 24(2): 1799. https://doi.org/10.3390/ijms24021799
  • Vermaas JV, Hardy DJ, Stone JE, Tajkhorshid E, Kohlmeyer A (2016) TopoGromacs: Automated Topology Conversion from CHARMM to GROMACS within VMD. Journal of Chemical Information and Modeling 56(6): 1112–1116. https://doi.org/10.1021/acs.jcim.6b00103
  • WHO (2020) Global Tuberculosis Report 2020. World Health Organisation.
  • Widawati (2022) Efektivitas Ekstrak Kulit dan Batang Lannea coromandelica sebagai Bahan Pengawet Anti Jamur Schizophyllum commune FRIES di bawah bimbingan Ira Taskirawati dan Syahidah.
  • Zhou Q-L, Tan Z-H, Wang H-X, Chen D-J, Ke X-R, Zhu Z-X, Wang H-F (2021) The complete plastome of a cultivar of Lannea coromandelica. Mitochondrial DNA, Part B, Resources 6: 3386–3387. https://doi.org/10.1080/23802359.2021.1998803
login to comment