Preprint
Article

In Silico Drug Re-Purposing Studies for the Discovery of Novel MbtA Inhibitors

This version is not peer-reviewed.

Submitted:

27 September 2023

Posted:

28 September 2023

You are already at the latest version

A peer-reviewed article of this preprint also exists.

Abstract
Tuberculosis (TB) continues to pose a global health challenge, exacerbated by the rise of drug-resistant strains. The development of new TB therapies is an arduous and time-consuming process. To expedite the discovery of effective treatments, computational structure-based drug re-purposing has emerged as a promising strategy. From this perspective, conditionally essential targets present a valuable opportunity, and the mycobactin biosynthesis pathway stands out as a prime example highlighting the intricate response of Mycobacterium tuberculosis (Mtb) to changes in iron availability. This study focuses on the re-purposing and revival of FDA-approved drugs (library) as potential inhibitors of MbtA, a crucial enzyme in mycobactin biosynthesis in Mtb conserved among all species of mycobacteria. Literature suggests this pathway to be associated with drug efflux pumps, which potentially contribute to drug resistance. This makes it a potential target for antitubercular drug discovery. Herein we utilized cheminformatics and structure-based drug repurposing approaches, viz., molecular docking, dynamics, and PCA analysis, to decode the intermolecular interactions and binding affinity of the FDA-reported molecules against MbtA. The virtual screening revealed ten molecules with significant binding affinities and interactions with MbtA. These drugs, originally designed for different therapeutic indications (4: antiviral, 3: anticancer, 1: CYP450 inhibitor, 1: ACE inhibitor, and 1: leukotriene antagonist), are repurposed as potential MbtA inhibitors. Furthermore, our study explores the binding modes and interactions between these drugs and MbtA, shedding light on the structural basis of their inhibitory potential. Principal component analysis highlighted significant motions in MbtA-bound ligands, emphasizing the stability of the top Protein-Ligand Complexes (PLCs). This computational approach provides a swift and cost-effective method for identifying new MbtA inhibitors, which can subsequently undergo validation through experimental assays. This streamlined process is facilitated by the fact that these compounds are already FDA-approved and have established safety and efficacy profiles. This study has the potential to lay the groundwork for addressing the urgent global health challenge at hand, specifically in the context of combating Antimicrobial Resistance (AMR) and Tuberculosis (TB).
Keywords: 
Subject: 
Medicine and Pharmacology  -   Pharmacology and Toxicology

1. Introduction

For centuries, Tuberculosis (TB) has haunted humanity as an infectious disease caused by the lethal pathogen Mycobacterium tuberculosis (Mtb). Today, this bacterium has escalated into a global pandemic that wreaks havoc worldwide[1]. The World Health Organization (WHO) released a 2021 report stating that TB caused 1.6 million deaths (including 187000 people with HIV) and affected more than 10.6 million people worldwide in 2021, signifying its dire impact on the mortality and morbidity rates globally[2]. As per the WHO's latest estimates, one in four individuals currently carry an established TB infection. Mtb is a severe bacterial infection that primarily affects the lungs but can also affect other parts of the body. Its ability to transition between respiring and non-respiring conditions without losing viability as a result of various enzymatic reactions leads to challenges in tackling TB infection[3]. It is a global health issue that requires improved living conditions, better access to healthcare, and increased awareness and prevention efforts to combat its spread. United Nations has been actively involved in efforts to combat tuberculosis (TB) globally. The UN's Sustainable Development Goals include a target to end the TB epidemic by 2030[2].
The escalation in global infections has led to drug discovery and development becoming central areas of focus within the scientific community. Antitubercular drug discovery began with the discovery of p-amino salicylic acid (1943), Streptomycin (1944), Isoniazid and Pyrazinamide (1952), Ethambutol (1961), Rifampicin (1963). Following this, there was a forty-year research pause primarily due to a lack of worldwide financing, resistant existing drug targets, unviable new drug targets, and clinical trial failures of novel medications. In an effort to mitigate the atrocities associated with infectious TB, a few clinically approved drugs were launched: Bedaquiline (2012), Delamanid (2014), and pretomanid (2019)[4]. These drugs continue to be the only choice of treatment for drug-resistant TB to date. Recently, TB has become an epidemic due to the surge in numerous resistance cases, particularly DR-TB, MDR-TB, XDR-TB, and TDR-TB cases[5,6]. These drug-resistant cases seem to be potentially incurable with available therapies. Moreover, the existence of co-infections pre- and post-COVID era has worsened the current scenario in tacking TB[7,8,9]. These factors necessitate the search for novel antitubercular agents acting via novel mechanisms of action on novel targets. This could open avenues for tackling drug resistance. Figure 1 outlines the timeline of drug discovery in TB alongside the present situation.
In this sense, CET (conditionally essential target) based drug design offers an excellent opportunity[10]. This encompasses focusing on a process that is conditionally necessary for mycobacteria, supporting their enzymatic functions and aiding in bacterial colonization, proliferation, and growth. Mycobacteria rely on the essential process of iron acquisition, vital for their biochemical machinery, from human host sources[11,12,13]. Although serum iron is scarce (10-24M), abundant iron is found in extracellular bound form. Mycobacteria have evolved various mechanisms to counter iron scarcity, notably synthesizing, secreting, and re-uptaking small molecules termed mycobactins/siderophores/iron chelators[14,15]. This siderophore machinery becomes active under iron-deficient conditions, a critical factor in tuberculosis pathogenesis[16]. The mycobactin biosynthesis pathway involves non-ribosomal peptide synthetase and polyketide synthase (NRPS-PKS) assembly chains[17]. Building blocks like salicylate and modified lysine are synthesized externally and linked to NRPS-PKS enzymes, as shown in Figure 3. Amino acids are added to growing mycobactin molecules. Mycobactins are hexadentate ligands with a tripodal structure, comprising an o-hydroxyphenyloxazoline system linked to an acylated hydroxylysine residue attached to a 3-hydroxybutyric acid unit, further tethered to a cyclized hydroxylysine, forming a seven-membered lactam ring[18]. Mtb maintains iron homeostasis, regulating uptake, utilization, and storage[19]. Under iron stress, genes for siderophore biosynthesis (mbt), export (MmpL4/5-MmpS4/5), and import (IrtAB) are upregulated[20]. Siderophore synthesis involves proteins like HupB, IdeR, and Lsr2. Chelating agents/mycobactins are released in the intracellular space. Carboxymycobactins are exported via inner-membrane transporters MmpL4/5 with MmpS4/5 adaptors. Mycobactin is exported via extracellular vesicles. Carboxymycobactins transfer chelated Fe3+ to mycobactin via HupB, acting as a receptor[21]. Once iron is chelated, it's internalized as a ferri-carboxymycobactin/ferri-mycobactin complex and imported independently via iron-regulated ABC transporter IrtAB[22]. IrtAB is a unique inner membrane heterodimer coupling Fe3+-carboxymycobactins import and iron assimilation through IrtA's cytosolic domain, functioning as a flavin reductase, converting Fe3+ to Fe2+[23]. Mtb stores iron in bacterioferritin (BfrA) and ferritin (BfrB) to maintain iron balance[24,25]. Iron toxicity is prevented by the type VII Esx-3 secretion system, although its precise role in iron uptake is unclear[26,27]. Carboxymycobactins and mycobactin, once deferrated, are recycled and exported via MmpL4/5. From this machinery, it is evident that the absence of the mycobactin gene inhibits growth in iron deficiency, highlighting the significance of the mbt gene cluster[28]. Figure 2 presents a highlight of the above-mentioned process.
Until now, researchers have explored merely four enzymes - MbtI, MbtA, MbtM, and PPTase - within the mycobactin biosynthesis pathway as potential drug targets against tuberculosis, and the advances have been published by our team in the form of a review [4]. Our primary focus lies on MbtA (salicyl-AMP ligase) due to its role in catalyzing the initial step of mycobactin biosynthesis and its critical role in bacterial growth and virulence. In the mycobactin megasynthase gene cluster, a crucial bimodular system comprising salicyl-AMP ligase (MbtA) and phenyloxazoline synthetase (MbtB) drives the production of the heterocyclic oxazoline segment within the mycobactin structure. MbtA plays a pivotal role in initiating mycobactin biosynthesis through two sequential reactions: first, by catalyzing the activation of salicylic acid, forming salicyl-adenosine monophosphate (Sal-AMP), and then by transferring Sal-AMP to the phosphopantetheinylation domain of MbtB. MbtB, an acyl carrier protein, acts as the adjoining enzyme to MbtA in the NRPS-PKS cluster. The ensuing MbtB-bound salicyl-thioester is condensed with ser-anchored to MbtB's C-terminal phosphopantetheinylation site - resulting in the formation of an amide. Through several subsequent steps, this amide is cyclized to generate the 2-hydroxyphenyloxazoline core structure of mycobactins. This MbtA-MbtB bimodular system is an indispensable and conserved element of the biochemical machinery employed by pathogenic mycobacterial species to survive in iron-deficient conditions within host macrophages, making it a promising target for therapeutic intervention[29]. Despite the existence of a few inhibitors (nucleoside and non-nucleoside) targeting MbtA[7,30,31], the quest to discover a more suitable inhibitor persists, and furthermore, the study of MbtA inhibitors is currently in its preliminary stages.
Figure 3. The biosynthetic pathway of mycobacterial siderophores, mycobactins & carboxymycobactins.
Figure 3. The biosynthetic pathway of mycobacterial siderophores, mycobactins & carboxymycobactins.
Preprints 86286 g003
We are pioneering drug repurposing for this target (MbtA) to aid global scientific efforts against tuberculosis and antimicrobial resistance. Our focus centers on repurposing FDA-approved drug compounds that target mycobacterial salicyl-AMP ligase (MbtA), a critical enzyme within the mycobactin biosynthetic pathway. Through the utilization of structure-based drug design and virtual screening, we systematically evaluated the entire molecular library of FDA-reported drugs to discern their binding and molecular interactions. This extensive assessment allowed us to gather valuable insights into the selective binding of FDA drugs to MbtA, enabling the assessment of their potential as therapeutic agents. To delve deeper, molecular dynamics simulations were employed to probe the molecular underpinnings of inhibition, exploring the binding modes and the stability of the complex formed between inhibitors and the target enzyme. These in silico analyses provide an initial and comprehensive understanding of the binding mechanisms and the potential of FDA molecules as MbtA inhibitors. This approach provides deeper insights into the individual amino acids' contributions to ligand binding and the presence of crucial protein-ligand interactions during simulations, along with frame analysis and PCA analysis. Such findings hold the promise of guiding the development of novel therapeutic agents that target MbtA, thereby paving the way for potential applications in combatting resistance challenges. The outcomes of this study hold the potential for informing the reformation and assessment of new potential medications to combat tuberculosis and antimicrobial resistance. Moreover, these findings contribute to our broader comprehension of inhibitory mechanisms. Insights gained from this study would increase our preparedness against TB, which can evolve into a deadly pandemic. Our team also has plans to carry out biological evaluations of the top-ranking hit compounds in order to assess their inhibitory activity against MbtA and as efflux pump modulators. These efforts would also contribute to addressing the global challenge of coinfections.

2. Materials and Methods

2.1. System Specifications and Software Employed

All molecular simulations were performed on a workstation with the following specifications: Operating System: Ubuntu 22.04 LTS, 64-bit, Processor: Intel® Core™ i5-12400 CPU@2.30 GHz processor, RAM: 16 GB, and Graphics: 8 GB Nvidia GeForce RTX 3050 GPU. The software employed for the study was: (i) The ligand database (FDA-Approved Drugs) for carrying out virtual screening was downloaded from the web server (https://zinc.docking.org/)[32,33]. The file named fda.smi was saved for further molecular simulation studies. (ii) The co-crystallized protein structure of MbtA was obtained from the Protein Data Bank (PDB) and Alpha Fold database[34,35]. (iii) Molecular docking and virtual screening were performed using AutoDock-GPU in the Google Colab web server (https://colab.research.google.com/). The best conformers of the ligands were generated using an in-house Python script, write_largest_cluster_ligand.py, located in Autodock 4.2.6[36]. The summary of the virtual screening was generated using an in-house Python script: summarize_docking.py. (iv) Protein-ligand complexes (PLCs) were generated using AutoDock 4.2.6[37]. To carry out the molecular dynamics simulations (MDS), the GROMACS (GROningen MAchine for Chemical Simulations) version 2022.4 was utilized[38]. (v) All 2D and 3D visualizations were done using PyMOL (molecular visualization system) and LigPlot+ (academic version)[39,40].

2.2. Ligand Preparation

To conduct the virtual screening by molecular docking, the fda.smi file containing the details of ligands downloaded from the ZINC database (https://zinc.docking.org/) was used[41]. Using Open Babel: The Open-Source Chemistry Toolbox v. 2.4.0, the ligand library was converted to individual ligand.pdb file format[42]. Then, the ligands (ligand.pdb) were imported into AutoDock-GPU in Google Colab, and their energy was minimized. This process aided in the determination of bond orders and the addition of hydrogens to the ligands required for the docking studies. The individual output files (ligand.pdbqt) containing the best conformations of the ligands were used for further docking-based screening studies.

2.3. Protein Preparation

The protein structure of Salicyl-AMP ligase/salicyl-S-ArCP synthetase (565 amino acids) was obtained from the Alpha Fold Protein Structure Database (https://alphafold.com/). The source organism was Mycobacterium tuberculosis (strain ATCC 25618 / H37Rv) with protein code P71716 (MBTA_MYCTU). The PDB-BLAST search outcomes demonstrated that Mtb-MbtA shares sequence identity with related proteins, such as the crystal structure of DhbE bound with 2,3-dihydroxybenzoic acid (DHB)-adenylate (PDB code 1MDB) from Bacillus subtilis and Mycobacterium smegmatis MbtA apo structure (PDB ID: 5KEI)[43,44]. The active site of the MbtA protein was identified by analogy with these homologous structures and was prepared using the protein tab in AutoDock 4.2.6 program by MGLTools 1.5.6 for docking. Preparation steps involved: (i) deletion of water molecules, (ii) addition of polar hydrogens, (iii) merging non-polar hydrogens, (iii) specifying AD4 atom type. And (iv) addition of gasteiger charges. Finally, the protein was energy minimized and saved in protein.pdbqt format for protein-ligand docking.

2.4. Identification of Binding Site and Receptor Grid Generation

Locating the active site is an essential parameter for protein-ligand docking. The probable active site in the model was identified by aligning it with the template of the ligand-bound DHB adenylating enzyme DhbE (PDB ID: 1MDB). This enabled a comparison of the active site residues in both structures, which was essential for conducting the following docking studies. The grid was generated based on the position of the co-crystallized ligand present in the active site of the selected protein (PDB: 1MDB). Subsequently, the MbtA protein was aligned with this ligand's position, and a grid box was formed around the centroid of the co-crystallized ligand. The receptor atoms' van der Waals radius was adjusted to 1.00 Å, along with a partial atomic charge of 0.25. This approach was chosen due to the conserved nature of the active site, given the homologous structure of both proteins. The input file (protein.gpf) was generated, and AutoGrid executable file was run to generate the output (protein.glg) file. The exact grid coordinates were further used in virtual screening, as presented in Table 1.

2.5. Docking-Based Virtual Screening Studies

Prior to commencing virtual screening, a validation of the docking program was performed. The co-crystallized ligand of 1MDB was re-docked into the active site of the MbtA protein (homology). Subsequently, the root mean square deviation (RMSD) was calculated. Then, docking-based screening of the FDA-reported molecules (fda.smi) was performed in the active site of the MbtA receptor (by analogy with the related structure PDB: 1MDB) to evaluate the binding modes and inhibitory profile in AutoDock-GPU in Google Colab[45]. Similar grid parameters were employed as in validation studies. The script is available as AutoDock_GPU_VJ.ipynb on GitHub (https://github.com/). The following docking parameters were established: a population size of 150 ensured accuracy within a reasonable computation time, a maximum of 2,500,000 energy evaluations were conducted to explore ample conformational space, and the genetic algorithms were capped at 27,000 generations for convergence. A total of 200 runs were executed to achieve a reliable sampling of conformational space. The visual representation (2D and 3D) of the docked structures was achieved using LigPlot+ v.2.2 and PyMOL.

2.6. Molecular Dynamics Studies

Utilizing the GROMACS package (version 2022.4, single precision, GPU enabled) along with the CHARMM force field, molecular dynamics simulations were performed[38]. The optimization of the MbtA protein structure involved the removal of water molecules, the addition of hydrogen atoms, and energy minimization through the steepest descent algorithm (5ns). TIP3P water molecules were employed to solvate the system within a cubic box, maintaining a protein-box edge distance of 10 Å. Neutralization of the system was achieved by introducing counterions (Na+ or Cl) utilizing the gmx_genion module. Simulations were executed within an NPT and NVT ensemble, adhering to periodic boundary conditions, operating at a temperature of 300 K and a pressure of 1 atm. Long-range electrostatic interactions were computed by the particle mesh Ewald method (cutoff distance: 12Å)[46]. The simulations were run for 300 ns with coordinates and velocities saved at 10ps intervals for subsequent analysis. The md300.xtc, md300.tpr, and md300.edr files generated were employed to generate RMSD, RMSF, ROG, and SASA data for analysis. Based on the stability trends observed in the RMSD graph of the protein-ligand complexes, a particular frame was selected for in-depth analysis. This frame serves as the foundation for a more comprehensive discussion on the essential residues involved in ligand binding with a focus on H-bonds.

2.7. Principal Component Analysis (PCA)

Principal Component Analysis (PCA) has been utilized to highlight the primary motions observed in protein trajectories (MD simulation) when bound to ligands. This analysis reveals that protein dynamics involve changes in molecular structure and conformation over time. PCA, a robust multivariate statistical technique, effectively reduces the dimensionality needed to describe protein dynamics through a systematic decomposition process, prioritizing observed motions from the largest to the smallest spatial scales[47]. The application of PCA to protein trajectories is referred to as Essential Dynamics (ED), aiming to extract the "essential" motions from sampled conformations. In practical terms, substantial large-scale motions can hinder the discernment of smaller-scale motions due to their significant amplitude in atomic displacements. It is noteworthy that these larger-scale motions often hold the utmost biological relevance in the realm of protein dynamics. Herein, PCA analysis was performed for all seven PLCs under study to explore molecular motion through Molecular Dynamics (MD) trajectories. To eliminate the translational and rotational motion of the molecule, a "least squares fit" to the reference structure was applied. This procedure involved obtaining a covariance matrix by linearly transforming Cartesian coordinate space. The matrix was then diagonalized to yield a set of eigenvectors that effectively represent the molecule's motion. Each eigenvector's associated eigenvalue indicated its energy contribution to the motion. Moreover, projecting the trajectory onto an eigenvector provided insights into the "time-dependent motions" exhibited by specific components within distinct vibrational modes. The temporal average of this projection helped discern the contribution of atomic vibrational components to this coordinated motion. To generate the eigenvectors and eigenvalues for the trajectory, the covariance matrix was computed and diagonalized using the integrated GROMACS utility "g_covar." Additionally, the "g_anaeig" tool was employed for examining and visualizing the eigenvectors[48]. In summary, for PCA analysis, a combination of techniques and tools, including least squares fitting, "g_covar," and "g_anaeig," integrated within GROMACS, were employed for a comprehensive analysis and visualization of molecular dynamics.

3. Results

3.1. Molecular Docking Studies on MbtA

3.1.1. Validation of Docking Procedure

The studies on the crystal structure of MbtA were validated through redocking, which resulted in a binding energy of -6.41 kcal/mol, a Ki of 20.06 µM, and a reference RMSD of 1.52 Å. Small fluctuations in RMSD (0-3Å) are acceptable for small globular macromolecules. Figure 4 presents the overlay conformations of the internal ligand with their co-crystallized conformation.

3.1.2. Virtual Screening of FDA-Reported Library Through Molecular Docking

This study utilizes the principles of structure-based drug design to identify potential molecules to tackle TB and AMR (antimicrobial resistance). Furthermore, we utilized the structure of MbtA (Salicyl-AMP ligase) to identify suitable binding agents with established safety profiles. For this purpose, FDA-approved drugs served as an ideal starting point, enabling the repurposing of safe and approved medications for combating tuberculosis. The virtual screening uncovered numerous molecules with favorable negative binding energies, including essential parameters such as docking scores, ligand efficiency, lipophilicity, hydrogen bonding interactions for each drug, and various other vital parameters. However, we opted to choose the top 10 hit molecules (results detailed in Table 2 and Table 3) with the highest negative binding energies. A higher negative binding energy indicates a stronger affinity for the active site. A bar graph representation of the same is shown in Figure 5 for a better understanding. Upon scrutinizing the molecules with the highest rankings, a notable observation emerged: among the top 10 molecules, four molecules (a_391: Saquinavir, a_85: Ritonavir, a_472: Lopinavir, and a_1276: Indinavir) were found to be protease inhibitors widely employed in antiviral therapies, notably recognized for their effectiveness against HIV and even could be used as boosters to antiviral therapies, 3 molecules possessing anti-cancer properties (a_821: Carfilzomib, a_1338: Venetoclax, and a_1388: Neratinib), one molecule as CYP450 inhibitor (a_617: Cobicistat), one ACE inhibitor (a_827: Candesartan), and one molecule as leukotriene antagonist (a_797: Zafirlukast (Accolate)). Thus, these molecules, stemming from various significant pharmacological classes, could potentially pave the way for exploring combination drug development aimed at tackling dual coinfections. To gain insights into the crucial interacting motifs, the ten compounds set was divided into groups on the basis of biological activity.
Interaction analysis of the top-scoring compounds obtained through virtual screening: Interpretation of molecular docking outcomes relies on specific descriptors, including binding energy, electrostatic energy, hydrogen bonding, van der Waals energy, and solvation energies[57]. Binding energy serves as a quantitative measure to compare and prioritize different ligands or potential drug candidates. It represents the overall energetic change associated with forming a stable complex between the ligand and receptor. Ligands with lower binding energies indicate a more favorable, strong, and specific interaction, suggesting a stronger binding affinity with the target receptor and vice-versa. Binding energies are often decomposed into various energy terms, as mentioned above, and analyzing these individual components provides a deeper understanding of the specific types of interactions driving the binding process. Electrostatic interactions (negative values) influence the optimal binding orientation of the ligand within the receptor's active site, leading to conformational changes in both the ligand and the receptor, thereby impacting the overall fit and stability of the complex. They also influence the solvent effects, charged residues, and specificity of binding, followed by structural rearrangements. Hydrogen bonds are formed when a hydrogen atom of the ligand's functional group interacts with electronegative atoms, such as oxygen or nitrogen, in the receptor (amino acids). Their presence, number, and strength contribute to the stability and specificity of the ligand-receptor complex, aiding in the accurate orientation of the ligand, and facilitating optimal interactions with the receptor's active site residues. Van der Waals energy quantifies the favorable interactions between hydrophobic portions of the ligand and receptor. These interactions contribute to the snug fitting of the ligand into the receptor's binding pocket, optimizing binding affinity and specificity. Solvation energies account for the energetics of solvating molecules, and they help simulate the ligand-receptor interactions in a more realistic physiological context, forming interactions with the receptor. They contribute to the overall binding free energy and provide insights into the balance between hydrophobic and hydrophilic interactions. To summarize, the assessment of these values during docking simulations provides a comprehensive understanding of the ligand's compatibility with the receptor's binding site. The detailed interactions of the top ten ligands have been presented in Table 4 and Figure 6a-j.
  • a_617: Cobicistat – MbtA complex (Figure 6a)
Cobicistat is a Cytochrome P450 3A inhibitor used as a pharmacokinetic enhancer in combination with certain HIV-1 protease inhibitors (PIs) to improve their effectiveness. This showed the highest binding energy of -16.69 kcal/mol, which signifies its effective binding in the active site of MbtA. Three H-bonds were made: Asn258, Thr462, and Arg451. (i) One from the OH4 of the carbonyl oxygen of a_617 to the nitrogen of Asn285 (OH4Carbonyl oxygen – NHAsn258 = 2.80Å), (ii) N7H of the thiazole ring of a_617 to the nitrogen of Thr462 (N7HThiazole Ring – NHThr462 = 3.06Å), and O1H of the carbonyl oxygen of a_617 to the oxygen (OG1) of Thr462 (O1HCarbonyl oxygen – OHThr462 = 2.98Å), and (iii) N4H of the thiazole ring of a_617 to the nitrogen2 (guanidine group) of Arg451 (N4HThiazole Ring – NH2Arg451 = 3.15Å), O2H of the carbonyl oxygen of a_617 to the nitrogen2 (guanidine group) of Arg451 (O2HCarbonyl oxygen – NH2Arg451 = 2.48Å), and O3 of the morpholine ring of a_617 to the nitrogen1 (guanidine group) of Arg451 (O3Morpholine ring – NH1Arg451 = 3.10Å).
2.
a_391: Saquinavir – MbtA complex (Figure 6b)
Saquinavir, an inhibitor of HIV protease, exhibited a binding energy of -16.33 kcal/mol. It established six H-bonds with active site residues: Asn258, Thr462, Arg451, Gly460, Ala356, & Phe259. (i) Carbonyl oxygen of N-tert-butylformamide (a_391) to the α-amino group (NH) of Asn258 (ON-tert-butylformamide – NHAsn258 = 2.87Å), (ii) hydroxyl oxygen (attached to octahydroisoquinolin) of a_391 to the oxygen of Thr462 (O2HOctahydroisoquinolin – OHThr462 = 2.81Å), (iii) N of the quinoline ring of a_391 to the nitrogen2 (guanidino group) of Arg451 (NQuinoline ring – NH2Arg451 = 3.21Å) and Oxygen of the oxopropyl-quinoline of a_617 to the nitrogen1 (guanidino group) of Arg451 (OOxopropyl-quinoline – NH1Arg451 = 3.23Å), (iv) hydroxyl oxygen (attached to octahydroisoquinolin) of a_391 to the Oxygen of Gly460 (O2HOctahydroisoquinolin – ONH2-CH2-COOH (Gly460) = 3.05Å), (v) Oxygen attached to quinoline-2-carboxamide (a_391) to the amino group (NH) of Ala356 (OQuinoline-2-carboxamide – NHAla356 = 2.92Å), and (vi) Carbonyl oxygen of N-tert-butylformamide (a_617) to the amino group (NH) of Phe259 (O1N-tert-butylformamide – NHPhe259 = 3.07Å).
3.
a_821: Carfilzomib – MbtA complex (Figure 6c)
Carfilzomib, acting as a proteasome inhibitor, exhibited a binding energy of -16.08 kcal/mol with three hydrogen bonds: Gly330, Thr462, & Arg451. (i) N1H of N-methylacetamide in a_821 to the OH group of Gly330 (NHN-methylacetamide – OHGly330 = 2.95Å), (ii) Oxygen of formamido-N-methylacetamide in a_821 to the side chain hydroxyl group of Thr462 (O5Formamido-N-methylacetamide – OHThr462 = 3.21Å), and (iii) Oxygen2 of N-methylacetamide of a_821 to the nitrogen2 (guanidino group) of Arg451 (O2N-methylacetamide – NH2Arg451 = 2.76Å)
4.
a_827: Candesartan – MbtA complex (Figure 6d)
Candesartan, an angiotensin II receptor antagonist, demonstrated a binding energy of -15.82 kcal/mol. It established five H-bonds with active site residues: Asn258, Arg451, Asp436, Gly354, & His257. (i) Oxygen of cyclohexyl hydrogen carbonate part in a_827 to the NH of side chain carboxamide of Asn258 (O4Cyclohexyl hydrogen carbonate – NHAsn258 = 2.96Å) and carbonyl oxygen (O2) in a_827 to the NH of α-amino group of Asn258 (O2 – NHAsn258 = 2.82Å), (ii) N6 of the tetrazole ring of a_827 to the nitrogen1 (guanidino group) of Arg451 (N6Tetrazole ring – NH1Arg451 = 3.06Å) and to the nitrogen2 (guanidino group) of Arg451 (N6Tetrazole ring – NH2Arg451 = 3.16Å), (iii) N5 of the tetrazole ring of a_827 to the acidic side chain (CH2COOH) of Asp436 (N5Tetrazole ring – OHAsp436 = 2.72Å), (iv) N1 of dihydro-1H-imidazole of a_827 to the NH of Gly354 (N1Dihydro-1H-imidazole – NHGly354 = 3.01Å), (v) Oxygen of cyclohexanol in a_827 to the NH of imidazole side chain of His257 (O6Cyclohexanol – NHHis257 = 2.87Å).
5.
a_85: Ritonavir – MbtA complex (Figure 6e)
Ritonavir, an HIV protease inhibitor, displayed a binding energy of -14.84 kcal/mol, engaging in four hydrogen bonds with residues Ala356, Val212, Gly354, and His257. (i) Carbonyl oxygen of N-phenethyl acetamide (a_85) to the amino group (NH) of Ala356 (ON-phenethyl acetamide – NHAla356 = 2.88Å), (ii) Oxygen5 in a_85 to the carbonyl oxygen of α-carboxylic acid group of Val212 (O5a_827– Oα-carboxylic acid:Val212 = 2.80Å), (iii) N2 of azanecarboxamide of a_85 to the carbonyl oxygen of Gly354 (N1Azanecarboxamide – OGly354 = 2.96Å), (iv) Oxygen5 in a_85 to the NH of imidazole side chain of His257 (O5a_827 – NHHis257 = 3.08Å).
6.
a_472: Lopinavir – MbtA complex (Figure 6f)
Lopinavir, an HIV protease inhibitor employed in the treatment of HIV infection, exhibited a binding energy of -14.92 kcal/mol and formed four hydrogen bonds with active site residues Gly354, Gly460, Thr462, and Arg451. (i) N1H of N-ethylpropionamide of a_617 to the carbonyl oxygen of Gly354 (N1HN-ethylpropionamide – OGly354 = 2.75Å), (ii) hydroxyl oxygen3 (attached to N-(2-hydroxypropyl)acetamide) of a_472 to the Oxygen of Gly460 (O3HN-(2-hydroxypropyl)acetamide – ONH2-CH2-COOH (Gly460) = 2.92Å), (iii) hydroxyl oxygen3 (attached to N-(2-hydroxypropyl)acetamide) of a_472 to the side chain hydroxyl group of Thr462 (O3HN-(2-hydroxypropyl)acetamide – OHThr462 = 3.21Å) and NH of tetrahydro-pyrimidin-2-one of a_472 to the hydroxyl of carboxyl group of Thr462 (NHTetrahydro-pyrimidin-2-one – OHThr462 = 3.16Å), and (iv) Carbonyl oxygen of N-methylacetamide of a_472 to the nitrogen2 (guanidino group) of Arg451 (ON-methylacetamide – NH2Arg451 = 2.86Å),
7.
a_1276: Indinavir – MbtA complex (Figure 6g)
Indinavir, an integral component of highly active antiretroviral therapy for treating HIV/AIDS, exhibited a strong binding affinity of -14.80 kcal/mol and displayed two hydrogen bonds with active site residues Gly354 and Arg451. (i) NH of piperazine-2-carboxamide of a_1276 to the carbonyl oxygen of Gly354 (NHPiperazine-2-carboxamide – OGly354 = 2.66Å) and (ii) Carbonyl oxygen of piperazine-2-carboxamide of a_1276 to the nitrogen2 (guanidino group) of Arg451 (OPiperazine-2-carboxamide – NH2Arg451 = 2.63Å) and nitrogen of pyridine of a_1276 to the nitrogen1 (guanidino group) of Arg451 (OPyridine – NH1Arg451 = 2.93Å).
8.
a_1338: Venetoclax – MbtA complex (Figure 6h)
Venetoclax, a B-cell lymphoma-2 (BCL-2) inhibitor known for its anti-apoptotic role, exhibited a robust binding to MbtA with a binding energy of -14.68 kcal/mol and established three hydrogen bonds with residues Val352, Thr462, and Phe259. (i) Nitrogen of benzenesulfonamide of a_1338 to the hydroxyl group of α-carboxylic acid group of Val352 (NBenzenesulfonamide – OHVal352 = 3.17Å), (ii) Oxygen of tetrahydro-2H-pyran of a_1338 to the side chain hydroxyl group of Thr462 (OTetrahydro-2H-pyran – OHThr462 = 2.83Å), and (iii) Carbonyl oxygen of nitrobenzene of a_1338 to the side chain amino group (NH) of Phe259 (O6Nitrobenzene – NHPhe259 = 3.07Å).
9.
a_797: Zafirlukast (Accolate) – MbtA complex (Figure 6i)
Zafirlukast, a leukotriene receptor antagonist, demonstrated binding to MbtA with a binding energy of -14.50 kcal/mol and formed three hydrogen bonds with active site residues Asn258, His257, and Arg451. (i) Carbonyl oxygen of cyclopentyl methylcarbamate (a_797) to the α-amino group (NH) of Asn258 (OCyclopentyl methylcarbamate – NHAsn258 = 3.07Å), (ii) Carbonyl oxygen of cyclopentyl methylcarbamate (a_797) to the NH of the imidazole side chain of His257 (OCyclopentyl methylcarbamate – NHHis257 = 3.04Å), and (iii) Carbonyl oxygen of benzenesulfonamide of a_797 to the NH (guanidino group) of Arg451 (O6Benzenesulfonamide - NHArg451 = 2.87Å) and hydroxyl of benzimidic acid of a_797 to the nitrogen2 (guanidino group) of Arg451 (OBenzimidic acid – NH2Arg451 = 2.88Å).
10.
a_1388: Neratinib – MbtA complex (Figure 6j)
Neratinib, recognized for its anticancer properties as a tyrosine kinase inhibitor, displayed a binding energy of -15.20 kcal/mol while establishing three hydrogen bonds with active site residues His257, Phe259, & Thr462. (i) Oxygen attached to dihydroquinoline of a_1388 to the NH of imidazole side chain of His257 (ODihydroquinoline – NHHis257 = 2.88Å), (ii) Carbonyl oxygen of N-methylacrylamide of a_1388 to the side chain amino group (NH) of Phe259 (O2N-methylacrylamide – NHPhe259 = 2.94Å), and (iii) Bridging Nitrogen (N3) of a_1388 to the side chain hydroxyl group of Thr462 (N3Linker bridge – OHThr462 = 2.94Å).
The active site residues (amino acids) Asn258, Thr462, Arg451, Gly460, Ala356, Phe259, Gly330, Asp436, Gly354, His257, Val212, and Val352 are crucially showing significant contribution towards binding of the ligands in the active site of MbtA (salicyl AMP-ligase) by forming hydrogen bonds with external ligands, particularly in the context of tuberculosis (TB) and mycobactin biosynthesis[58]. These interactions contribute to substrate recognition, binding, and catalysis within the active site. Asn258 is likely involved in forming hydrogen bonds with ligands, aiding in their precise orientation and stabilization within the active site. Thr462 and Arg451 form multiple hydrogen bonds to anchor ligands and stabilize their binding conformation within the active site. The presence of Gly460 can influence the conformation and flexibility of nearby residues, potentially affecting ligand interactions. Ala356 may contribute to the overall structural stability of the active site and provide a hydrophobic environment for ligand binding. Phe259 and His257 can form π-π interactions or other types of hydrogen bonds with ligands, aiding in their recognition, stabilization, and binding specificity. Gly330, being similar to Gly460, may influence the flexibility of nearby residues and contribute to the overall dynamics of ligand binding. Asp436 can form hydrogen bonds with ligands and assist in positioning them for catalysis or recognition. Gly354, as with other glycine residues, can impact the overall flexibility of the active site region, influencing ligand interactions. Val212 and Val352 potentially create a hydrophobic environment within the active site, aiding in ligand binding through hydrogen bonds and hydrophobic interactions while also assisting in the proper positioning of ligands. All other active site residues play a role in facilitating hydrophobic interactions by protein folding, contributing significantly to the stability and biological activity of proteins. They enable proteins to minimize their surface area, thereby reducing unfavorable interactions with water. Hence, in the context of MbtA, the formation of hydrogen bonds and hydrophobic residues between these amino acids and external ligands contributes to the specificity and efficiency of substrate/ligand binding and catalytic processes.
Table 4. The table details hydrogen bonding and hydrophobic interactions of the top ten ligands in MbtA's active site residues.
Table 4. The table details hydrogen bonding and hydrophobic interactions of the top ten ligands in MbtA's active site residues.
Sl. No. Code Docking interactions with active site amino acid residues H-bond distance (Å)
1 a_617 H-bond- Asn258, Thr462, & Arg451
Hydrophobic- Val448, Asp436, Lys332, Phe353, Gly450, Gly330, Gly354, Leu360, Val352, Phe259, Cys263, Leu253, Ser213, Glu461, Ala356, Gly460, His257, Gly214, Met355, Thr216, Glu357, & Cys457
2.80, [3.06, 2.98], & [3.15, 2.48, 3.10]
2 a_391 H-bond- Asn258, Thr462, Arg451, Gly460, Ala356, & Phe259
Hydrophobic- His257, Leu253, Glu461, Val302, Leu360, Val352, Gly354, Gly329, Gln376, Ser331, Gly330, Phe353, Lys332, Val448, Met355, Thr216, Glu357, Val212, Gly214, & Ser213
2.87, 2.81, [3.23, 3.21], 3.05, 2.92, & [2.87, 3.07]
3 a_821 H-bond- Gly330, Thr462, & Arg451
Hydrophobic- Ala254, His257, Glu461, His523, Leu253, Leu126, His129, Gly214, Ser213, Val212, Met355, Ala356, Ser434, Glu357, Tyr432, Thr216, Val352, Lys332, Ser331, Phe353, Gly354, Gly460, & Ala459
2.95, 3.21, & 2.76
4 a_827 H-bond- Asn258, Arg451, Asp436, Gly354, & His257
Hydrophobic-Val212, Thr462, Glu461, Ala356, Val352, Gly330, Phe259, Leu360, Ser331, Phe353, Ser434, Tyr432, Met355, Tyr415, Gly214, Glu357, Thr216, & Ser213
[2.82, 2.96], [3.16, 3.06], 2.72, 3.01, 2.87
5 a_85 H-bond- Ala356, Val212, Gly354, & His257
Hydrophobic- Gly214, Val352, Gly330, Met355, Thr462, Phe259, Asn258, Pro260, Leu253, Phe353, Val302, Ser331, Gly460, Glu357, Leu126, Glu461, Arg451, Ser213, Val448, Cys457, & Asp436
2.88, 2.80, 2.96, & 3.08
6 a_472 H-bond- Gly354, Gly460, Thr462, & Arg451
Hydrophobic- Asp436, Val448, Phe353, Met355, Lys332, Gly330, Ser331, Val352, Phe259, Leu360, Gly214, Val212, His129, Ser213, Leu126, Ala356, Glu461, Ala459, His257, & Leu253
2.75, 2.92, 3.16, & 2.86
7 a_1276 H-bond- Gly354 & Arg451
Hydrophobic- Asp436, Phe353, Leu360, Phe259, Gly329, Val352, Asn258, Glu461, Val302, Gly330, His257, Ala254, Thr462, Leu253, Gly460, Ala356, Glu357, Tyr432, Thr216, Gly214, & Ser434
2.66, & [2.63, 2.93]
8 a_1338 H-bond- Val352, Thr462, & Phe259
Hydrophobic- Gly460, Met355, Ser434, Tyr432, Arg451, Asp436, Cys457, Ile456, Gly450, Val455, Glu493, Glu334, Ser331, Val448, Lys332, Gly330, Phe353, Gly354, Gly329, Asn258, His257, Ala356, & Leu360
3.17, 2.83, & 3.15
9 a_797 H-bond- Asn258, His257, & Arg451
Hydrophobic- Glu461, Gly214, Gly460, Ser213, Thr462, Val212, Lys332, Asp436, Ser331, Cys457, Gly450, Val448, Phe353, Gly354, Val352, Gly330, Cys263, Leu360, Phe259, & Ala356
3.07, 3.04, & [2.87, 2.88]
10 a_1388 H-bond- His257, Phe259, & Thr462
Hydrophobic- Arg451, Val352, Phe353, Gly329, Ser331, Val448, Gln376, Lys332, Gly330, Asp436, Met355, Gly354, Ser213, Gly214, Glu357, Ala356, Gly256, Pro260, Gly460, Val212, Glu461, Val302, Leu253, & Asn258
2.88, 2.94, & 2.94

3.2. Molecular Dynamics Simulation Study

To gain a more profound insight into the binding mechanism and validate the ligand-protein complexes, we conducted 300 ns molecular dynamics (MD) simulations to analyze the biophysical interactions between the docked complexes (seven in number). These selections encompassed a wide range of compound classes of interest, including protease inhibitors utilized in HIV therapy, CYP inhibitors, ACE inhibitors, and anticancer agents. These identified compounds have the potential to act as dual inhibitors, exhibiting both binding affinity to and inhibition of MbtA. The evaluation included Root Mean Square Deviation (RMSD), Root Mean Square Fluctuation (RMSF), and the percentage of interactions between the ligands and protein atoms. The retention of interactions observed during docking was also verified over the 300 ns MD simulation. Conformational changes were measured using RMSD values, where the initial frame served as the reference backbone, and the deviation at the end was compared to the starting point. The system achieved stabilization and equilibration when RMSD values ranged from 0.1 to 1.5 Å. Lower complex stability was associated with more significant fluctuations. This analysis allowed us to ascertain the reliability of the ligand-protein interactions and provided insights into the dynamic behavior of the complexes during the simulation.
RMSD analysis: In complex a_617-MbtA, initially, the RMSD was volatile till 100ns, but then it was stable throughout the simulation with an RMSD of 0.6nm (protein) and 0.4nm (ligand). For complex a_391-MbtA, the complex was stable throughout the simulation with a constant RMSD of 0.45nm (protein) and 0.35nm (ligand). There were mild fluctuations in the initial frames (0-30ns) but eventually, the complex stabilized around fixed values. For complex a_821-MbtA, the complex was stable between 50-200ns at a fixed RMSD of 0.5nm. Fluctuations were observed initially in both protein and ligand, but after 200ns, the ligand displayed fluctuations with RMSD in the range of 0.7-1.1nm. For complex a_827-MbtA, the protein fluctuated significantly for the first 110ns. However, it then attained stability and had an RMSD of 0.5nm. The ligand was stable throughout the course of the simulation. However, it displayed stability at two RMSD values, 0-160ns at 0.25nm and 160-300ns at 0.6nm. Overall, the complex was stable from 160-300ns at an RMSD of 0.6nm. For complex a_85-MbtA, the ligand was stable throughout the simulation (0-300ns) at a constant RMSD of 0.25nm while the protein fluctuated a lot initially (0-160ns: 0.3-0.55nm) and towards the end of the simulation (280-300ns: 0.45-0.55nm) and eventually in mid-segment it gained stability (160-280ns) at a RMSD of 0.4nm. In total, the complex was stable from 160-280ns. For complex a_472-MbtA, it was very surprising to observe that the protein was stable throughout the course of the simulation (0-300ns) at a constant RMSD of 0.5nm. The ligand was stable initially (0-90ns) at an RMSD of 0.5, but eventually, it had major fluctuations throughout the rest of the simulation (90-300ns). Overall, the complex was stable from 0-90ns at an RMSD of 0.5, and for complex a_1276-MbtA, there were initial fluctuations in RMSD for both protein and ligand, but eventually, from 50-160ns both the protein and ligand had a constant RMSD of 0.6nm post which (160-300ns) there were minor fluctuations but the RMSD window was large (protein:0.45-0.8 and ligand: 0.4-0.8nm). The identified trends indicate sustained stability across all protein-ligand complexes during the simulation, punctuated by minor fluctuations, as demonstrated in Figure 7. Among these, PLC (a), (b), (c), and (g) displayed remarkable stability throughout the simulation duration. Meanwhile, PLC (d), (e), and (f) exhibited intermittent fluctuations; nevertheless, they maintained a commendable level of overall stability.
RMSF analysis: RMSF is a measure of the fluctuation of atomic positions in a biomolecular system over the course of a simulation and is determined by the root mean square deviation of each atom's position from its average position over the course of a simulation. The peaks visible in Figure 8 correspond to areas of the protein that exhibited the highest levels of fluctuation throughout the simulation. Based on our experience, the N- and C-terminal tails of the protein tend to be the most dynamic and variable regions. Distinct variations were notable in the RMSF plots for C-alpha residues ranging from 60-70, 180-220, 280-320, and 450-550 across all seven complexes. These observed variations are likely indicative of the spatial adjustments occurring in the active site residues upon binding of the ligands. In the case of complex a_827-MbtA, noticeable fluctuations were evident, as depicted by a distinct peak in the 60-70 residue range due to their flexibility. For complex a_472-MbtA, a substantial number of variations were observed in residues 450-550. Nevertheless, it is worth noting that all complexes exhibited satisfactory RMSF values, indicating their overall stability. It's also evident that specific active site residues exhibit a degree of flexibility in the presence of ligand-bound complexes. This suggests that the formation of complexes with inhibitors contributes to the stabilization of the active site.
Rg analysis: The radius of gyration (Rg) plot is a graphical representation that illustrates the compactness or size of a molecule during simulation, as shown in Figure 9. The radius of gyration measures the average distance of each atom in a molecule from its center of mass. For a protein, a lower radius of gyration indicates a more compact and folded structure, while a higher value suggests a more extended or unfolded structure. In our investigation, all protein-ligand complexes, with the exception of the a_1276-MbtA complex, exhibited an initial expansion in the size of the protein-ligand assembly. This phenomenon pointed to conformational changes occurring during the simulation, which were subsequently followed by a reduction and stabilization phase around 50 ns. Beyond this point, intermittent spikes were evident, suggesting a certain degree of structural stability and the emergence of enduring interactions between the protein and the ligand. However, the a_1276-MbtA complex displayed commendable stability initially, maintaining this state until approximately 60 ns. Following this period, there was a sharp rise in spikes and an increase in the size of the protein-ligand assembly after 150 ns, indicative of a substantial conformational alteration. Despite this, the complex managed to attain stability again, albeit at a higher value. The oscillations seen in the Rg plot signify instances of unfolding or alterations in the conformation, while the plateaus indicate periods of steady structural arrangements. Notably, the protein-ligand complex initiates the simulation with observable conformational changes, as evidenced by the initial expansion in size. Subsequently, the complex experiences stabilization, evident through a decrease and steadying of the Rg values, interspersed with episodes of fluctuations. This pattern strongly suggests the development of enduring interactions between the protein and the ligand, progressively solidifying as the simulation unfolds.
SASA Analysis: SASA quantifies the surface area of a molecule that is accessible to solvent molecules, providing valuable insights into the molecule's exposure to its surrounding environment. In our study, apart from the a_1276-MbtA complex, all other complexes exhibited a comparable exposure pattern. They demonstrated an initial rise between 0-25 ns, followed by a consistently stable graph with minimal spikes and fluctuations throughout the entire simulation period. The general upward trend in SASA across a_1276-MbtA complex complexes suggests that the protein-ligand interactions gradually expose themselves to the solvent over the course of the simulation, while this was constant for other ligands as they converged around fixed values. The rise in SASA could potentially signify unfolding or conformational changes, whereas a decrease might imply the adoption of more compact and structurally stable conformations. This observation attests to the stability of all the protein-ligand complexes over the simulation time, as shown in Figure 10.
Time Frame Analysis with a focus on H-bond and active site residues: Frame analysis in molecular dynamics refers to the systematic examination and evaluation of the molecular structure and behavior at specific time intervals or frames during a simulation. It plays a crucial role in understanding the dynamics, stability, and interactions of molecules over time. This total simulation period is divided into discrete timeframes or frames, where each frame represents a specific point in time. Upon analysis of the highest-ranking molecules from Table 3, time frame analysis was conducted by extracting a snapshot of the frame from the region of most excellent stability, as indicated by the RMSD (Root Mean Square Deviation) graph. The gmx_trjconv module was employed to take a snapshot. It is evident that, during the simulation, this highly stable frame represents the interactions that are most likely to remain constant. These stable interactions are crucial for the effective binding of the ligand to the receptor. The detailed analysis is as follows: (a) a_617-MbtA complex made one H-bond: Nitrogen of thiazole of a_617 to the α-amino group (NH) of Gly214 (N6Thiazole – NHGly214 = 3.07Å), (b) a_391-MbtA complex made four H-bonds: (i) Carbonyl oxygen attached to quinoline-2-carboxamide to the α-amino group (NH) of Gly214 (O2Quinoline-2-carboxamide – NHGly215 = 3.02Å), (ii) Carbonyl oxygen attached to quinoline-2-carboxamide to the α-amino group (NH) of Gly215 (O2Quinoline-2-carboxamide – NHGly215 = 3.16Å), and (iii) Carbonyl oxygen (attached to octahydroisoquinolin) of a_391 to the NH of Thr462 (O1Octahydroisoquinolin – NHThr462 = 3.07Å) and N3H of a_391 to the carbonyl oxygen of Thr462 (N3Ha_391 – OThr462 = 3.02Å); this interaction was conserved as in docking, (iv) Carbonyl oxygen (attached to octahydroisoquinolin) of a_391 to the OH of carboxylic acid of Glu461 (O1Octahydroisoquinolin – OHGlu461 = 2.81Å) and N1H of a_391 to the OH of carboxylic acid of Glu461 (N1Ha_391 – OHGlu461 = 2.86Å), and (c) a_821-MbtA complex made two H-bonds: (i) NH of N-methylacetamide in a_821 to the OH group of Gly330 (NHN-methylacetamide – OHGly330 = 2.79Å), and (ii) NH of N-methyl-2-morpholinoacetamide in a_821 to the OH group of Gly460 (NHN-methyl-2-morpholinoacetamide – OHGly460 = 2.93Å). As an observation, it was seen in the top three molecules, the majority of the hydrogen bonds were established with glycine residue. The presence of Glycine can influence the conformation and flexibility of nearby residues, potentially affecting ligand interactions due to its small and structurally flexible nature. It also facilitates specific interactions with ligands. However, the specific role of glycine residues in the active site of MbtA can vary depending on their positions within the protein's primary structure. For (d) a_827-MbtA complex, there were two H-bonds: (i) N4 of the tetrazole ring of a_827 to the acidic side chain (CH2COOH) of Asp436 (N4Tetrazole ring – OHAsp436 = 2.53Å), and (ii) Nitrogen of imidazole in a_827 to the NH of imidazole side chain of His257 (NImidazole – NHHis257 = 2.84Å). The a_85-MbtA complex displayed two H-bonds; (i) Carbonyl oxygen of N-phenethyl acetamide (a_85) to the amino group (NH) of Ala356 (ON-phenethyl acetamide – NHAla356 = 2.66Å), and (ii) Carbonyl oxygen of a_85 to the amino group (NH) of Phe259 (O4a_391 – NHPhe259 = 3.22Å). Complex a_472-MbtA made one H-bond; NH of tetrahydro-pyrimidin-2-one of a_472 to the hydroxyl of the carboxyl group of Thr462 (NHTetrahydro-pyrimidin-2-one – OHThr462 = 2.86Å), and a_1276-MbtA complex showed no H-bonding. The majority of interactions remained consistent across all these four protein-ligand complexes. Figure 11 displays the interactions made by the ligands in complex with MbtA for the most stable frame during MD simulation, as evidenced by the RMSD graph. The research pertaining to the exploration of the crystal structure of MbtA as an antitubercular target is at a very initial stage and calls for much more evidence-based literature. In summary, the frame analysis of our top-scoring PLCs provides a systematic examination of molecular properties and interactions at discrete time intervals for the most stable frame.

4. Principal Component Analysis (PCA)

Principal Component Analysis (PCA) investigated the primary modes of motion within the simulation trajectory. To comprehend the primary structural variations elucidated by each MbtA-ligand complex, we generated PCA graphs. The collective motion represented by the first two principal components (PCs) and 2D projections of PC1 and PC2 were plotted for all the Protein-Ligand Complexes (PLCs), as depicted in Figure 12. As an observation, it was seen that the complex a_85-MbtA, a_821-MbtA, and a_1276-MbtA express the compact clusters in the conformational spaces that range from – 10 to 10; for the latter two PLCs, there were slight scatterings across the configurational space (Figure 12 e,c,g). In the MD trajectory of a_617-MbtA and a_391-MbtA, complex PC1 and PC2 (top two modes) showed a varied constant distribution across the configurational space with little variation in the conformational space and was widely grouped in the range of – 15 to 5. (Figure 12a,b). In the case of a_827-MbtA and a_472-MbtA complex, the graph displayed a varied distribution across the configurational space with slightly more variation as evident from the widely grouped data in the range of -5 to 10 (Figure 12d,f).
In these plots (Figure 13), distinct color transitions were evident at various points, signifying shifts between different conformations induced by ligand binding. Each data point represents an individual frame within the trajectory. Throughout the 300 ns simulation, the systems generally maintained compactness, with minor deviations observed for all PLCs. These fluctuations could be attributed to temporary disruptions in hydrogen bonds or the accommodation of ligands within the active site. These analyses underscore the considerable flexibility of the protein structure when interacting with the seven hit compounds during the initial simulation phases, with a gradual reduction in flexibility as the simulation progressed. Furthermore, the contribution percentage in eigenmodes consistently decreased, indicating a tendency towards local structural stabilization in each complex. These motions primarily stemmed from the interactions of the docked compounds within the active site of the protein, reinforcing the formation of robust complexes for each protein-ligand pair.
Root Mean Square Fluctuation (RMSF) nm of eigenvectors-backbone in PCA of molecular dynamics trajectories refers to the measure of atomic fluctuations or deviations from their average positions within a protein over a period of time during a simulation, thereby providing insights into the flexibility and mobility of specific regions or residues within the structure. Herein, as an observation, it was seen that vector 1 displayed fluctuations in the atom number range 200-250, and all other vectors displayed significant fluctuations in the atom range 1250-1750 except vector3 for all PLCs. This information provides insights into the relative stability of different parts of the molecule throughout the simulation. It also assists in understanding structural dynamics and may offer valuable clues about functional or binding sites. The minimal fluctuations observed indicate that all seven PLCs maintain stability within the active site of MbtA, suggesting their potential contribution to effective inhibition. The graphical representation of RMSF has been presented in Figure 14.
In summary, PCA analysis revealed that the interconnected motion within the MbtA conformation reflects both the rigidity and significant variations induced at the active site as a consequence of ligand binding throughout the simulation. These analyses strongly indicate the stability of the top seven selected compounds within the MbtA active site, which, in turn, hinders the essential motions required for enzymatic activity. Multiple factors, including RMSD, RMSF, Rg, protein-ligand interactions, PCA studies, and numerous molecular docking scores, substantiate this inhibition.

5. Discussion

Tackling TB is crucial for improving public health, ensuring global health security, and advancing toward a world where infectious diseases are less of a threat to human well-being. The mycobactin biosynthesis machinery is conserved among all the species of mycobacteria and presents a novel target for antitubercular drug development. Our primary focus lies on MbtA (salicyl-AMP ligase) due to its pivotal role in catalyzing the initial step of mycobactin biosynthesis and its critical role in bacterial growth and virulence. Moreover, the modulation of efflux pumps could also help tackle the emergence of AMR in TB. Hence, to rediscover a drug with alternative therapeutic applications as a MbtA inhibitor, we conducted virtual screening using the FDA-approved drug library against the active pocket of the MbtA receptor. Our virtual screening yielded promising results, indicating strong binding affinities. Based on the calculated binding energies, we selected the top ten molecules for further investigation. Surprisingly, the top hits encompassed a diverse range of significant pharmacological classes: four molecules (a_391: Saquinavir, a_85: Ritonavir, a_472: Lopinavir, and a_1276: Indinavir) are well-known protease inhibitors widely utilized in antiviral therapies, three molecules are anti-cancer agents (a_821: Carfilzomib, a_1338: Venetoclax, and a_1388: Neratinib), one functions as a CYP450 inhibitor (a_617: Cobicistat), one as an ACE inhibitor (a_827: Candesartan), and one as a leukotriene antagonist (a_797: Zafirlukast, or Accolate). These compounds exhibited docking scores ranging from -16.69 to -14.14 kcal/mol, indicative of strong and specific interactions with the target receptor. The interaction between the ligands and the amino acid residues within the active site of MbtA plays a critical role in substrate recognition, binding, and catalysis. Key active site amino acid residues, including Asn258, Thr462, Arg451, Gly460, Ala356, Phe259, Gly330, Asp436, Gly354, His257, Val212, and Val352, participate in hydrogen bonding and various other interactions with the ligands. These interactions anchor the ligands and stabilize their binding conformation within the active site, potentially influencing the conformation and flexibility of neighbouring residues. To gain a deeper understanding of the binding characteristics, we conducted molecular dynamics simulations (300ns) for the top seven molecules. Analysis of the Root Mean Square Deviation (RMSD) revealed the stability of all seven Protein-Ligand Complexes (PLCs). The protein exhibited RMSD values within the range of 0.2-0.6 nm, while the ligands displayed RMSD values ranging from 0.2 to 1.5 nm, suggesting stable and consistent interactions within the complexes. The stability of the PLCs was confirmed through the comparison of protein and ligands RMSD graphs. These graphs exhibited consistent, smooth curves with minimal spikes for all PLCs, indicating their stable behaviour throughout the simulation (average 100ns). The analysis of Root Mean Square Fluctuation (RMSF) showed minimal atomic position fluctuations within the PLCs. However, distinct variations were observed in the RMSF plots for C-alpha residues, particularly in the ranges of 60-70, 180-220, 280-320, and 450-550, across all seven complexes. These variations likely signify spatial adjustments occurring within the active site residues upon ligand binding. Furthermore, it is evident that specific active site residues display a degree of flexibility when in the presence of ligand-bound complexes. This suggests that the formation of these complexes with inhibitors plays a role in stabilizing the active site. The analysis of the Radius of Gyration (Rg) for all PLCs revealed a notable pattern in the size of the protein-ligand assembly during the simulation. Initially, there was an expansion in the size, indicating conformational changes taking place. This phase was followed by a reduction and stabilization around the 50 ns mark. Beyond this point, intermittent spikes were observed, signifying a degree of structural stability and the emergence of enduring interactions between the protein and the ligand. The analysis of Solvent Accessible Surface Area (SASA) revealed a consistent exposure pattern among all Protein-Ligand Complexes (PLCs). Initially, there was an upward trend observed between 0-25 ns, which was succeeded by a persistently stable graph with minimal spikes and fluctuations throughout the entire simulation duration. The increase in SASA values during the initial phase may suggest unfolding or conformational alterations within the complexes. Conversely, a subsequent decrease could indicate the adoption of more compact and structurally stable conformations. This observation provides evidence of the stability maintained by all the protein-ligand complexes throughout the simulation period. To gain more insight into the binding modes and to identify structural changes in molecules or complexes, a particular time frame analysis (trajectory snapshot) was carried out for each PLCs from the most stable region of the RMSD plot. Results revealed that intermolecular interactions change over time, including atomic positions, velocities, energies, bond lengths, angles, dihedral angles, and solvent interactions, among others. This also highlights the system's dynamics, stability, and equilibrium behavior (H-bond interaction) over time. The stability of the Protein-Ligand Complexes (PLCs) could be attributed to the establishment of hydrogen bonds with specific residues. Notably, residues Gly214, Gly215, Glu461, Thr462, Gly460, Gly330, Asp436, His257, Phe259, Ala356, and Thr216 are involved in forming these hydrogen bonds, which likely play a significant role in ensuring stable binding mechanisms. Furthermore, Principal Component Analysis (PCA) substantiated the earlier findings, highlighting the noteworthy motions and stabilities within the MbtA protein when bound to the top seven ligands. In summary, these seven molecules show promise as potential MbtA inhibitors due to their strong affinity for the active site and their potential to act as substrates for enzyme inhibition. Hence this current study is solely based on the computational approaches to assess the inhibitory interaction of FDA-reported molecules against MbtA.

6. Conclusion

The mycobacterial machinery constitutes an intricately complex process, demanding immediate attention for drug discovery. The urgency stems from the emergence of diverse resistant strains and the waning efficacy of current drugs. Taking into account the severity of this bacterial infestation along with antimicrobial resistance, there is a need to design, discover, and rediscover drugs acting on a novel target via a novel mechanism of action having broad-spectrum antibacterial activity. Our study focuses on targeting the critical dual enzyme system, MbtA-MbtB, which plays a vital role in the biosynthesis of iron-chelating mycobacterial siderophores. Inhibiting these enzymes with a drug would not only disrupt the siderophore production machinery but also inhibit the operation of efflux pumps. Employing the concept of Drug repurposing, also known as drug repositioning or reprofiling, we evaluated the FDA library of drugs against the MbtA protein. We, fortunately, ended up rediscovering seven hit molecules that were initially developed for different therapeutic indications as MbtA inhibitors. We employed the concepts of virtual screening and also made a discussion on the binding modes and interactions between these drugs and MbtA, shedding light on the structural basis of their inhibitory potential. It is noteworthy that these molecules exhibit well-established ADMET profiles, further bolstering their suitability for repurposing. This innovative utilization of drug repurposing not only underscores the adaptability of existing pharmaceuticals but also offers a practical avenue for effectively addressing the complexities of co-infection scenarios in today's healthcare landscape.

7. Future Scope

As previously mentioned, it was quite remarkable to discover that among the top 10 molecules, four of them are protease inhibitors widely used in antiviral therapies, particularly effective against HIV. This suggests the potential for these compounds to act as dual inhibitors in tackling cases of co-infection. The future direction of our research involves evaluating these molecules for their inhibitory activity against MbtA in both GAST and GAST-Fe Media against M. smegmatis and M. tuberculosis. Additionally, we plan to conduct a Universal CAS siderophore inhibition assay in M. smegmatis. These assays represent a crucial step towards uncovering the dual inhibitory nature of these FDA-reported molecules. Having established safety and efficacy profiles, these drugs have the potential to become remarkable assets in the near future for addressing dual coinfections and combatting antimicrobial resistance. We are currently in the process of procuring these Active Pharmaceutical Ingredients (APIs) and will soon commence the studies. We anticipate publishing the next version of this research in the near future.

Supplementary Materials

Not applicable.

Author Contributions

“Conceptualization, VJ and GR; methodology, GR; software, GR and AB; validation, GR; formal analysis, GR and VJ; investigation, GR; resources, GR; data curation, GR; writing—original draft preparation, GR; writing—review and editing, GR and VJ; visualization, GR; supervision, VJ; project administration, VJ. All authors have read and agreed to the published version of the manuscript.”.

Funding

This research received no external funding.

Institutional Review Board Statement

Not applicable.

Informed Consent Statement

Not applicable.

Acknowledgments

The authors express their gratitude to the Department of Pharmaceutical Sciences and Technology at BIT Mesra for providing access to the computational facilities at the CADD laboratory and for granting the use of a workstation to conduct all molecular simulations. Gourav Rakshit is thankful to Birla Institute of Technology, Mesra, for providing financial assistance in the form of an Institute Research Fellowship (March 2021-2025). We extend our appreciation to all the open-source Software that has enabled us to conduct our studies with precision.

Conflicts of Interest

The authors declare no conflict of interest.

References

  1. Borah, P.; Deb, P.K.; Venugopala, K.N.; Al-Shar’i, N.A.; Singh, V.; Deka, S.; Srivastava, A.; Tiwari, V.; Mailavaram, R.P. Tuberculosis: An Update on Pathophysiology, Molecular Mechanisms of Drug Resistance, Newer Anti-TB Drugs, Treatment Regimens and Host-Directed Therapies. Curr Top Med Chem 2021, 21, 547–570. [Google Scholar] [CrossRef] [PubMed]
  2. World Health Organization Tuberculosis: Key Facts. Available online: https://www.who.int/news-room/fact-sheets/detail/tuberculosis (accessed on 18 September 2023).
  3. Peddireddy, V.; Doddam, S.N.; Ahmed, N. Mycobacterial Dormancy Systems and Host Responses in Tuberculosis. Front Immunol 2017, 8, 84. [Google Scholar] [CrossRef] [PubMed]
  4. Shyam, M.; Shilkar, D.; Verma, H.; Dev, A.; Sinha, B.N.; Brucoli, F.; Bhakta, S.; Jayaprakash, V. The Mycobactin Biosynthesis Pathway: A Prospective Therapeutic Target in the Battle against Tuberculosis. J Med Chem 2020, 64, 71–100. [Google Scholar] [CrossRef]
  5. Prasad, R.; Singh, A.; Balasubramanian, V.; Gupta, N. Extensively Drug-Resistant Tuberculosis in India: Current Evidence on Diagnosis & Management. Indian J Med Res 2017, 145, 271. [Google Scholar]
  6. Velayati, A.A.; Farnia, P.; Hoffner, S. Drug-Resistant Mycobacterium Tuberculosis: Epidemiology and Role of Morphological Alterations. J Glob Antimicrob Resist 2018, 12, 192–196. [Google Scholar] [CrossRef] [PubMed]
  7. Rakshit, G.; Jayaprakash, V. Tuberculosis and HIV Responses Threatened by NCOVID-19: A Situation Prompting an in Silico Investigation of Reported MbtA Inhibitors for Combined Inhibition of SARS-CoV-2 and HIV-TB Co-Infection. Struct Chem 2023, 34, 655–679. [Google Scholar] [CrossRef]
  8. Sibi, D.; Sethi Das, C.; Jibin, V.G.; Silvanose, C.D. Emerging Antibiotic Resistance in Post-COVID-19 Co-Infections. J Clin Med Case Reports 2023, 8, 8. [Google Scholar]
  9. Pawlowski, A.; Jansson, M.; Sköld, M.; Rottenberg, M.E.; Källenius, G. Tuberculosis and HIV Co-Infection. PLoS Pathog 2012, 8, e1002464. [Google Scholar] [CrossRef]
  10. Shyam, M.; Shilkar, D.; Rakshit, G.; Jayaprakash, V. Approaches for Targeting the Mycobactin Biosynthesis Pathway for Novel Anti-Tubercular Drug Discovery: Where We Stand. Expert Opin Drug Discov 2022, 17, 699–715. [Google Scholar] [CrossRef] [PubMed]
  11. Miethke, M.; Marahiel, M.A. Siderophore-Based Iron Acquisition and Pathogen Control. Microbiology and molecular biology reviews 2007, 71, 413–451. [Google Scholar] [CrossRef]
  12. Quadri, L.E.N. Iron Uptake in Mycobacteria. The mycobacterial cell envelope 2008, 167–184. [Google Scholar]
  13. Hammer, N.D.; Skaar, E.P. Molecular Mechanisms of Staphylococcus Aureus Iron Acquisition. Annu Rev Microbiol 2011, 65, 129–147. [Google Scholar] [CrossRef] [PubMed]
  14. Barclay, R.; Ewing, D.F.; Ratledge, C. Isolation, Identification, and Structural Analysis of the Mycobactins of Mycobacterium Avium, Mycobacterium Intracellulare, Mycobacterium Scrofulaceum, and Mycobacterium Paratuberculosis. J Bacteriol 1985, 164, 896–903. [Google Scholar] [CrossRef]
  15. Hall, R.M.; Ratledge, C. A Simple Method for the Production of Mycobactin, the Lipid-Soluble Siderophore, from Mycobacteria. FEMS Microbiol Lett 1982, 15, 133–136. [Google Scholar] [CrossRef]
  16. De Voss, J.J.; Rutter, K.; Schroeder, B.G.; Su, H.; Zhu, Y.; Barry III, C.E. The Salicylate-Derived Mycobactin Siderophores of Mycobacterium Tuberculosis Are Essential for Growth in Macrophages. Proceedings of the National Academy of Sciences 2000, 97, 1252–1257. [Google Scholar] [CrossRef]
  17. Barry, S.M.; Challis, G.L. Recent Advances in Siderophore Biosynthesis. Curr Opin Chem Biol 2009, 13, 205–215. [Google Scholar] [CrossRef]
  18. Snow, G. Mycobactins: Iron-Chelating Growth Factors from Mycobacteria. Bacteriol Rev 1970, 34, 99–125. [Google Scholar] [CrossRef]
  19. Sritharan, M. Iron Homeostasis in Mycobacterium Tuberculosis: Mechanistic Insights into Siderophore-Mediated Iron Uptake. J Bacteriol 2016, 198, 2399–2409. [Google Scholar] [CrossRef]
  20. Wells, R.M.; Jones, C.M.; Xi, Z.; Speer, A.; Danilchanka, O.; Doornbos, K.S.; Sun, P.; Wu, F.; Tian, C.; Niederweis, M. Discovery of a Siderophore Export System Essential for Virulence of Mycobacterium Tuberculosis. PLoS Pathog 2013, 9, e1003120. [Google Scholar] [CrossRef]
  21. Choudhury, M.; Koduru, T.N.; Kumar, N.; Salimi, S.; Desai, K.; Prabhu, N.P.; Sritharan, M. Iron Uptake and Transport by the Carboxymycobactin-Mycobactin Siderophore Machinery of Mycobacterium Tuberculosis Is Dependent on the Iron-Regulated Protein HupB. BioMetals 2021, 34, 511–528. [Google Scholar] [CrossRef]
  22. Rodriguez, G.M.; Smith, I. Identification of an ABC Transporter Required for Iron Acquisition and Virulence in Mycobacterium Tuberculosis. J Bacteriol 2006, 188, 424–430. [Google Scholar] [CrossRef] [PubMed]
  23. Ryndak, M.B.; Wang, S.; Smith, I.; Rodriguez, G.M. The Mycobacterium Tuberculosis High-Affinity Iron Importer, IrtA, Contains an FAD-Binding Domain. J Bacteriol 2010, 192, 861–869. [Google Scholar] [CrossRef] [PubMed]
  24. Gupta, V.; Gupta, R.K.; Khare, G.; Salunke, D.M.; Tyagi, A.K. Crystal Structure of Bfr A from Mycobacterium Tuberculosis: Incorporation of Selenomethionine Results in Cleavage and Demetallation of Haem. PLoS ONE 2009, 4, e8028. [Google Scholar] [CrossRef]
  25. Yao, H.; Wang, Y.; Lovell, S.; Kumar, R.; Ruvinsky, A.M.; Battaile, K.P.; Vakser, I.A.; Rivera, M. The Structure of the BfrB–Bfd Complex Reveals Protein–Protein Interactions Enabling Iron Release from Bacterioferritin. J Am Chem Soc 2012, 134, 13470–13481. [Google Scholar] [CrossRef] [PubMed]
  26. Siegrist, M.S.; Unnikrishnan, M.; McConnell, M.J.; Borowsky, M.; Cheng, T.-Y.; Siddiqi, N.; Fortune, S.M.; Moody, D.B.; Rubin, E.J. Mycobacterial Esx-3 Is Required for Mycobactin-Mediated Iron Acquisition. Proceedings of the National Academy of Sciences 2009, 106, 18792–18797. [Google Scholar] [CrossRef] [PubMed]
  27. Serafini, A.; Pisu, D.; Palù, G.; Rodriguez, G.M.; Manganelli, R. The ESX-3 Secretion System Is Necessary for Iron and Zinc Homeostasis in Mycobacterium Tuberculosis. PLoS ONE 2013, 8, e78351. [Google Scholar] [CrossRef]
  28. Jones, C.M.; Wells, R.M.; Madduri, A.V.R.; Renfrow, M.B.; Ratledge, C.; Moody, D.B.; Niederweis, M. Self-Poisoning of Mycobacterium Tuberculosis by Interrupting Siderophore Recycling. Proceedings of the National Academy of Sciences 2014, 111, 1945–1950. [Google Scholar] [CrossRef]
  29. Shyam, M.; Verma, H.; Bhattacharje, G.; Mukherjee, P.; Singh, S.; Kamilya, S.; Jalani, P.; Das, S.; Dasgupta, A.; Mondal, A. Mycobactin Analogues with Excellent Pharmacokinetic Profile Demonstrate Potent Antitubercular Specific Activity and Exceptional Efflux Pump Inhibition. J Med Chem 2022, 65, 234–256. [Google Scholar] [CrossRef] [PubMed]
  30. Stirrett, K.L.; Ferreras, J.A.; Jayaprakash, V.; Sinha, B.N.; Ren, T.; Quadri, L.E.N. Small Molecules with Structural Similarities to Siderophores as Novel Antimicrobials against Mycobacterium Tuberculosis and Yersinia Pestis. Bioorg Med Chem Lett 2008, 18, 2662–2668. [Google Scholar] [CrossRef]
  31. Ferreras, J.A.; Gupta, A.; Amin, N.D.; Basu, A.; Sinha, B.N.; Worgall, S.; Jayaprakash, V.; Quadri, L.E.N. Chemical Scaffolds with Structural Similarities to Siderophores of Nonribosomal Peptide-Polyketide Origin as Novel Antimicrobials against Mycobacterium Tuberculosis and Yersinia Pestis. Bioorg Med Chem Lett 2011, 21, 6533–6537. [Google Scholar] [CrossRef]
  32. Sterling, T.; Irwin, J.J. ZINC 15–Ligand Discovery for Everyone. J Chem Inf Model 2015, 55, 2324–2337. [Google Scholar] [CrossRef]
  33. Irwin, J.J.; Tang, K.G.; Young, J.; Dandarchuluun, C.; Wong, B.R.; Khurelbaatar, M.; Moroz, Y.S.; Mayfield, J.; Sayle, R.A. ZINC20—A Free Ultralarge-Scale Chemical Database for Ligand Discovery. J Chem Inf Model 2020, 60, 6065–6073. [Google Scholar] [CrossRef] [PubMed]
  34. Varadi, M.; Anyango, S.; Deshpande, M.; Nair, S.; Natassia, C.; Yordanova, G.; Yuan, D.; Stroe, O.; Wood, G.; Laydon, A.; et al. AlphaFold Protein Structure Database: Massively Expanding the Structural Coverage of Protein-Sequence Space with High-Accuracy Models. Nucleic Acids Res 2022, 50, D439–D444. [Google Scholar] [CrossRef] [PubMed]
  35. Berman, H.M.; Westbrook, J.; Feng, Z.; Gilliland, G.; Bhat, T.N.; Weissig, H.; Shindyalov, I.N.; Bourne, P.E. The Protein Data Bank. Nucleic Acids Res 2000, 28, 235–242. [Google Scholar] [CrossRef]
  36. Morris, G.M.; Goodsell, D.S.; Huey, R.; Hart, W.E.; Halliday, S.; Belew, R.; Olson, A.J. AutoDock. Automated docking of flexible ligands to receptor-User Guide, 2001. [Google Scholar]
  37. Rizvi, S.M.D.; Shakil, S.; Haneef, M. A Simple Click by Click Protocol to Perform Docking: Autodock 4.2 Made Easy for Non-Bioinformaticians. EXCLI J 2013, 12, 830–857. [Google Scholar]
  38. Van Der Spoel, D.; Lindahl, E.; Hess, B.; Groenhof, G.; Mark, A.E.; Berendsen, H.J.C. GROMACS: Fast, Flexible, and Free. J Comput Chem 2005, 26, 1701–1718. [Google Scholar] [CrossRef]
  39. DeLano, W.L. Pymol: An Open-Source Molecular Graphics Tool. CCP4 Newsl. Protein Crystallogr 2002, 40, 82–92. [Google Scholar]
  40. Laskowski, R.A.; Swindells, M.B. LigPlot+: Multiple Ligand–Protein Interaction Diagrams for Drug Discovery 2011.
  41. Irwin, J.J.; Shoichet, B.K. ZINC− a Free Database of Commercially Available Compounds for Virtual Screening. J Chem Inf Model 2005, 45, 177–182. [Google Scholar] [CrossRef]
  42. O’Boyle, N.M.; Banck, M.; James, C.A.; Morley, C.; Vandermeersch, T.; Hutchison, G.R. Open Babel: An Open Chemical Toolbox. J Cheminform 2011, 3, 1–14. [Google Scholar] [CrossRef]
  43. May, J.J.; Kessler, N.; Marahiel, M.A.; Stubbs, M.T. Crystal Structure of DhbE, an Archetype for Aryl Acid Activating Domains of Modular Nonribosomal Peptide Synthetases. Proceedings of the National Academy of Sciences 2002, 99, 12120–12125. [Google Scholar] [CrossRef]
  44. Vergnolle, O.; Xu, H.; Tufariello, J.M.; Favrot, L.; Malek, A.A.; Jacobs, W.R.; Blanchard, J.S. Post-Translational Acetylation of MbtA Modulates Mycobacterial Siderophore Biosynthesis. Journal of Biological Chemistry 2016, 291, 22315–22326. [Google Scholar] [CrossRef] [PubMed]
  45. Schieffer, G.; Peng, I. Accelerating Drug Discovery in AutoDock-GPU with Tensor Cores. In Proceedings of the European Conference on Parallel Processing; Springer; 2023; pp. 608–622. [Google Scholar]
  46. Kalibaeva, G.; Ferrario, M.; Ciccotti, G. Constant Pressure-Constant Temperature Molecular Dynamics: A Correct Constrained NPT Ensemble Using the Molecular Virial. Mol Phys 2003, 101, 765–778. [Google Scholar] [CrossRef]
  47. David, C.C.; Jacobs, D.J. Principal Component Analysis: A Method for Determining the Essential Dynamics of Proteins. Protein dynamics: Methods and protocols 2014, 193–226. [Google Scholar]
  48. Amadei, A.; Linssen, A.B.M.; Berendsen, H.J.C. Essential Dynamics of Proteins. Proteins: Structure, Function, and Bioinformatics 1993, 17, 412–425. [Google Scholar] [CrossRef]
  49. Deeks, E.D. Cobicistat: A Review of Its Use as a Pharmacokinetic Enhancer of Atazanavir and Darunavir in Patients with HIV-1 Infection. Drugs 2014, 74, 195–206. [Google Scholar] [CrossRef] [PubMed]
  50. Noble, S.; Faulds, D. Saquinavir: A Review of Its Pharmacology and Clinical Potential in the Management of HIV Infection. Drugs 1996, 52, 93–112. [Google Scholar] [CrossRef] [PubMed]
  51. Khan, M.L.; Stewart, A.K. Carfilzomib: A Novel Second-Generation Proteasome Inhibitor. Future oncology 2011, 7, 607–612. [Google Scholar] [CrossRef]
  52. Husain, A.; Azim, M.S.; Mitra, M.; Bhasin, P.S. A Review on Candesartan: Pharmacological and Pharmaceutical Profile. J Appl Pharm Sci 2011, 12–17. [Google Scholar]
  53. Mangum, E.M.; Graham, K.K. Lopinavir-ritonavir: A New Protease Inhibitor. Pharmacotherapy: The Journal of Human Pharmacology and Drug Therapy 2001, 21, 1352–1363. [Google Scholar] [CrossRef]
  54. Plosker, G.L.; Noble, S. Indinavir: A Review of Its Use in the Management of HIV Infection. Drugs 1999, 58, 1165–1203. [Google Scholar] [CrossRef]
  55. Chung, K.F.; Barnes, P.J. Zafirlukast (Accolate). Drugs Today (Barc) 1998, 34, 375–388. [Google Scholar] [CrossRef] [PubMed]
  56. Tiwari, S.R.; Mishra, P.; Abraham, J. Neratinib, a Novel HER2-Targeted Tyrosine Kinase Inhibitor. Clin Breast Cancer 2016, 16, 344–348. [Google Scholar] [CrossRef] [PubMed]
  57. Morris, G.M.; Lim-Wilby, M. Molecular Docking. Molecular modeling of proteins 2008, 365–382. [Google Scholar]
  58. Labello, N.P.; Bennett, E.M.; Ferguson, D.M.; Aldrich, C.C. Quantitative Three Dimensional Structure Linear Interaction Energy Model of 5′-O-[N-(Salicyl) Sulfamoyl] Adenosine and the Aryl Acid Adenylating Enzyme MbtA. J Med Chem 2008, 51, 7154–7160. [Google Scholar] [CrossRef] [PubMed]
Figure 1. The timeline of drug discovery in TB alongside the present situation.
Figure 1. The timeline of drug discovery in TB alongside the present situation.
Preprints 86286 g001
Figure 2. Mtb's Iron Response: when iron is scarce, Mtb boosts genes for siderophore synthesis.
Figure 2. Mtb's Iron Response: when iron is scarce, Mtb boosts genes for siderophore synthesis.
Preprints 86286 g002
Figure 4. The superimposed overlay conformation of the docked internal ligand (green) concerning its crystallized conformation (pink) obtained from the co-crystallized complex structure.
Figure 4. The superimposed overlay conformation of the docked internal ligand (green) concerning its crystallized conformation (pink) obtained from the co-crystallized complex structure.
Preprints 86286 g004
Figure 5. Graph comparing the lowest binding energies/docking scores and the number of conformations in the cluster (largest cluster) for top-scoring docked molecule in MbtA protein's active site.
Figure 5. Graph comparing the lowest binding energies/docking scores and the number of conformations in the cluster (largest cluster) for top-scoring docked molecule in MbtA protein's active site.
Preprints 86286 g005
Figure 6. The docked confirmation of (a) a_617, (b) a_391, (c) a_821, (d) a_827, (e) a_85, (f) a_472, (g) a_1276, (h) a_1338, (i) a_797, and (j) a_1388 in the active site of MbtA highlighting various interactions. In the plot, hydrophobic interactions are depicted as maroon spiked arcs, while hydrogen bonding interactions are shown as dashed green lines, with their lengths indicated in Å. The color scheme distinguishes various atoms and bonds: carbon atoms are black, oxygen atoms are red, sulphur as yellow, and nitrogen atoms are blue. Amino acid residue bonds are represented in brown, and ligands are colored violet.
Figure 6. The docked confirmation of (a) a_617, (b) a_391, (c) a_821, (d) a_827, (e) a_85, (f) a_472, (g) a_1276, (h) a_1338, (i) a_797, and (j) a_1388 in the active site of MbtA highlighting various interactions. In the plot, hydrophobic interactions are depicted as maroon spiked arcs, while hydrogen bonding interactions are shown as dashed green lines, with their lengths indicated in Å. The color scheme distinguishes various atoms and bonds: carbon atoms are black, oxygen atoms are red, sulphur as yellow, and nitrogen atoms are blue. Amino acid residue bonds are represented in brown, and ligands are colored violet.
Preprints 86286 g006
Figure 7. Root mean square deviation (RMSD) plot of MbtA protein (maroon) with ligand (a) a_617: Cobicistat, (b) a_391: Saquinavir, (c) a_821: Carfilzomib, (d) a_827: Candesartan, (e) a_85: Ritonavir, (f) a_472: Lopinavir, and (g) a_1276: Indinavir.
Figure 7. Root mean square deviation (RMSD) plot of MbtA protein (maroon) with ligand (a) a_617: Cobicistat, (b) a_391: Saquinavir, (c) a_821: Carfilzomib, (d) a_827: Candesartan, (e) a_85: Ritonavir, (f) a_472: Lopinavir, and (g) a_1276: Indinavir.
Preprints 86286 g007
Figure 8. Root mean square fluctuation (RMSF) plot of MbtA protein (green) with ligand (a) a_617: Cobicistat (blue), (b) a_391: Saquinavir (orange), (c) a_821: Carfilzomib (grey), (d) a_827: Candesartan (yellow), (e) a_85: Ritonavir (blue), (f) a_472: Lopinavir (green), and (g) a_1276: Indinavir (purple).
Figure 8. Root mean square fluctuation (RMSF) plot of MbtA protein (green) with ligand (a) a_617: Cobicistat (blue), (b) a_391: Saquinavir (orange), (c) a_821: Carfilzomib (grey), (d) a_827: Candesartan (yellow), (e) a_85: Ritonavir (blue), (f) a_472: Lopinavir (green), and (g) a_1276: Indinavir (purple).
Preprints 86286 g008
Figure 9. Radius of gyration (Rg) plot of MbtA protein (purple) with ligand (a) a_617: Cobicistat (blue), (b) a_391: Saquinavir (red), (c) a_821: Carfilzomib (yellow), (d) a_827: Candesartan (brown), (e) a_85: Ritonavir (black), (f) a_472: Lopinavir (green), and (g) a_1276: Indinavir (grey).
Figure 9. Radius of gyration (Rg) plot of MbtA protein (purple) with ligand (a) a_617: Cobicistat (blue), (b) a_391: Saquinavir (red), (c) a_821: Carfilzomib (yellow), (d) a_827: Candesartan (brown), (e) a_85: Ritonavir (black), (f) a_472: Lopinavir (green), and (g) a_1276: Indinavir (grey).
Preprints 86286 g009
Figure 10. Solvent accessible surface area (SASA) plot of MbtA protein (green) with ligand (a) a_617: Cobicistat (blue), (b) a_391: Saquinavir (red), (c) a_821: Carfilzomib (yellow), (d) a_827: Candesartan (brown), (e) a_85: Ritonavir (black), (f) a_472: Lopinavir (green), and (g) a_1276: Indinavir (grey).
Figure 10. Solvent accessible surface area (SASA) plot of MbtA protein (green) with ligand (a) a_617: Cobicistat (blue), (b) a_391: Saquinavir (red), (c) a_821: Carfilzomib (yellow), (d) a_827: Candesartan (brown), (e) a_85: Ritonavir (black), (f) a_472: Lopinavir (green), and (g) a_1276: Indinavir (grey).
Preprints 86286 g010
Figure 11. Interactions made by the ligands (a) a_617: Cobicistat, (b) a_391: Saquinavir, (c) a_821: Carfilzomib, (d) a_827: Candesartan, (e) a_85: Ritonavir, (f) a_472: Lopinavir, and (g) a_1276: Indinavir in complex with MbtA for the most stable frame during MD simulation as evidenced from RMSD graph. In the plot, hydrophobic interactions are depicted as maroon spiked arcs, while hydrogen bonding interactions are shown as dashed green lines, with their lengths indicated in Å. The color scheme distinguishes various atoms and bonds: carbon atoms are black, oxygen atoms are red, sulfur as yellow, and nitrogen atoms are blue. Amino acid residue bonds are represented in brown, and ligands are colored violet.
Figure 11. Interactions made by the ligands (a) a_617: Cobicistat, (b) a_391: Saquinavir, (c) a_821: Carfilzomib, (d) a_827: Candesartan, (e) a_85: Ritonavir, (f) a_472: Lopinavir, and (g) a_1276: Indinavir in complex with MbtA for the most stable frame during MD simulation as evidenced from RMSD graph. In the plot, hydrophobic interactions are depicted as maroon spiked arcs, while hydrogen bonding interactions are shown as dashed green lines, with their lengths indicated in Å. The color scheme distinguishes various atoms and bonds: carbon atoms are black, oxygen atoms are red, sulfur as yellow, and nitrogen atoms are blue. Amino acid residue bonds are represented in brown, and ligands are colored violet.
Preprints 86286 g011
Figure 12. Principal component analysis of protein-ligand complexes: the collective motion of (a) a_617: Cobicistat (yellow), (b) a_391: Saquinavir (green), (c) a_821: Carfilzomib (brown), (d) a_827: Candesartan (grey), (e) a_85: Ritonavir (red), (f) a_472: Lopinavir (blue), and (g) a_1276: Indinavir (black) with MbtA using projections of MD trajectories on two eigenvectors corresponding to the first two principal components.
Figure 12. Principal component analysis of protein-ligand complexes: the collective motion of (a) a_617: Cobicistat (yellow), (b) a_391: Saquinavir (green), (c) a_821: Carfilzomib (brown), (d) a_827: Candesartan (grey), (e) a_85: Ritonavir (red), (f) a_472: Lopinavir (blue), and (g) a_1276: Indinavir (black) with MbtA using projections of MD trajectories on two eigenvectors corresponding to the first two principal components.
Preprints 86286 g012
Figure 13. Principal component analysis of individual protein–ligand complexes: the collective motion of (a) a_617: Cobicistat, (b) a_391: Saquinavir, (c) a_821: Carfilzomib, (d) a_827: Candesartan, (e) a_85: Ritonavir, (f) a_472: Lopinavir, and (g) a_1276: Indinavir with MbtA using projections of MD trajectories on five various eigenvectors (1_2dproj: black, 2_2dproj: red, 3_2dproj: red, 4_2dproj: blue).
Figure 13. Principal component analysis of individual protein–ligand complexes: the collective motion of (a) a_617: Cobicistat, (b) a_391: Saquinavir, (c) a_821: Carfilzomib, (d) a_827: Candesartan, (e) a_85: Ritonavir, (f) a_472: Lopinavir, and (g) a_1276: Indinavir with MbtA using projections of MD trajectories on five various eigenvectors (1_2dproj: black, 2_2dproj: red, 3_2dproj: red, 4_2dproj: blue).
Preprints 86286 g013
Figure 14. Root mean square fluctuation (RMSF) nm of eigenvectors-backbone highlighting the stability of the protein in the presence of ligand (a) a_617: Cobicistat (blue), (b) a_391: Saquinavir (red), (c) a_821: Carfilzomib (yellow), (d) a_827: Candesartan (brown), (e) a_85: Ritonavir (black), (f) a_472: Lopinavir (green), and (g) a_1276: Indinavir (grey).
Figure 14. Root mean square fluctuation (RMSF) nm of eigenvectors-backbone highlighting the stability of the protein in the presence of ligand (a) a_617: Cobicistat (blue), (b) a_391: Saquinavir (red), (c) a_821: Carfilzomib (yellow), (d) a_827: Candesartan (brown), (e) a_85: Ritonavir (black), (f) a_472: Lopinavir (green), and (g) a_1276: Indinavir (grey).
Preprints 86286 g014
Table 1. Specifics of the grid parameter data utilized.
Table 1. Specifics of the grid parameter data utilized.
AutoDock 4.2.6
Protein Center grid box
dimensions (Å)
Spacing (Å) Coordinates for the center of the grid box
X-axis Y-axis Z-axis 0.375 X-axis Y-axis Z-axis
MbtA 50 50 50 -2.439 19.515 32.895
Table 2. Comprehensive tabular overview of the top ten hits resulting from docking-based virtual screening of FDA library against MbtA receptor.
Table 2. Comprehensive tabular overview of the top ten hits resulting from docking-based virtual screening of FDA library against MbtA receptor.
S. No Code Structure Compounds Details MOA and Therapeutic application
1 a_617 Preprints 86286 i001 ZINC ID: ZINC000085537014
Cobicistat
Formula: C40H53N7O5S2
Molar mass: 776.03 g·mol−1
Cobicistat (Tybost™) functions as a mechanism-based inhibitor of cytochrome P450 (CYP) 3A enzymes. It is approved in the EU as a pharmacokinetic enhancer, specifically a booster, of HIV-1 protease inhibitors (PIs) such as atazanavir and darunavir in adult patients[49]
2 a_391 Preprints 86286 i002 ZINC ID: ZINC000003914596
Saquinavir
Formula: C38H50N6O5
Molar mass: 670.855 g·mol−1
Saquinavir was the first protease inhibitor developed for HIV therapy. It acts by cleaving viral polypeptide chains into functional proteins during the late stages of HIV replication[50]
3 a_821 Preprints 86286 i003 ZINC ID: ZINC000049841054
Carfilzomib
Formula: C38H50N6O5
Molar mass: 670.855 g·mol−1
Carfilzomib (formerly PR-171) is a novel epoxyketone-based irreversible proteasome inhibitor. It specifically targets the chymotrypsin-like activity of the proteasome, preventing its normal function and causing a buildup of proteins within the cancer cells[51].
4 a_827 Preprints 86286 i004 ZINC ID: ZINC000004074875
Candesartan
Formula: C24H20N6O3
Molar mass: 440.463 g·mol−1
Candesartan blocks the vasoconstrictor and aldosterone-secreting effects of angiotensin II by selectively blocking the binding of angiotensin II to the AT1 receptor in many tissues, such as vascular smooth muscle and the adrenal gland[52].
5 a_85 Preprints 86286 i005 ZINC ID: ZINC000003944422
Ritonavir
Formula: C37H48N6O5S2
Molar mass: 720.95 g·mol−1
Ritonavir is a protease inhibitor and is used as booster with other protease inhibitors in the treatment of HIV infection. It prevents the proper maturation of the virus[53].
6 a_472 Preprints 86286 i006 ZINC ID: ZINC000003951740
Lopinavir
Formula: C37H48N4O5
Molar mass: 628.814 g·mol−1
Lopinavir is an antiretroviral medication primarily used in the treatment of HIV infection. It acts as a protease inhibitor in the final stages of the viral replication process thereby preventing the proper maturation of the virus leading to immature and non-functional viral particles[53].
7 a_1276 Preprints 86286 i007 ZINC ID: ZINC000022448696
Indinavir
Formula: C36H47N5O4
Molar mass: 613.803 g·mol−1
Indinavir is a protease inhibitor used as a component of highly active antiretroviral therapy to treat HIV/AIDS. It interferes with viral maturation, and contributes to the control of viral replication. It is also used as booster with other protease inhibitors[54].
8 a_1338 Preprints 86286 i008 ZINC ID: ZINC000150338755
Venetoclax
Formula: C45H50ClN7O7S
Molar mass: 867.320 g·mol−1
Venetoclax is B-cell lymphoma-2 (BCL-2) inhibitors. It binds to the BCL-2 protein, blocking its anti-apoptotic function
9 a_797 Preprints 86286 i009 ZINC ID: ZINC000000896717
Accolate
Formula: C31H33N3O6S
Molar mass: 575.210 g·mol−1
Zafirlukast (Accolate) is a leukotriene receptor antagonist. It inhibits the activity of leukotrienes by binding to their receptors found on various cells in the airways causing bronchodilation[55].
10 a_1388 Preprints 86286 i010 ZINC ID: ZINC000003916214
Neratinib
Formula: C30H29ClN6O3
Molar mass: 556.200 g·mol−1
Neratinib is a tyrosine kinase inhibitor used in the treatment of certain types of breast cancers (targeted therapy). It works by irreversibly inhibiting multiple receptor tyrosine kinases, including HER2, epidermal growth factor receptor (EGFR), and HER4[56].
Table 3. Detailed tabular representation of top 10 hits obtained through docking-based virtual screening from FDA-library.
Table 3. Detailed tabular representation of top 10 hits obtained through docking-based virtual screening from FDA-library.
S. No Ligand Runs Binding Energy (kcal/mol) ESTAT HB VDW DSOLV
1 a_617 200 -16.69 -0.9438 -2.5163 -23.9398 5.5597
2 a_391 200 -16.33 -0.6982 -2.18 -21.8025 4.9488
3 a_821 200 -16.08 -1.1546 -1.562 -23.9542 5.2432
4 a_827 200 -15.82 -0.8058 -1.9775 -21.7725 4.9785
5 a_85 200 -14.84 -0.3628 -2.1603 -20.6416 4.0052
6 a_472 200 -14.82 -0.8465 -2.6141 -20.8497 5.0782
7 a_1276 200 -14.80 -0.9278 -2.6608 -20.6792 5.1776
8 a_1338 200 -14.68 -1.3963 -1.6507 -27.4664 6.791
9 a_797 200 -14.50 -0.0758 -1.9036 -22.3757 4.7207
10 a_1388 200 -14.14 -0.0544 -1.7261 -21.6482 4.6123
Disclaimer/Publisher’s Note: The statements, opinions and data contained in all publications are solely those of the individual author(s) and contributor(s) and not of MDPI and/or the editor(s). MDPI and/or the editor(s) disclaim responsibility for any injury to people or property resulting from any ideas, methods, instructions or products referred to in the content.
Copyright: This open access article is published under a Creative Commons CC BY 4.0 license, which permit the free download, distribution, and reuse, provided that the author and preprint are cited in any reuse.
Alerts
Prerpints.org logo

Preprints.org is a free preprint server supported by MDPI in Basel, Switzerland.

Subscribe

© 2025 MDPI (Basel, Switzerland) unless otherwise stated