Open access peer-reviewed chapter

Molecular Informatics of Trypanothione Reductase of Leishmania major Reveals Novel Chromen-2-One Analogues as Potential Leishmanicides

Written By

Samuel K. Kwofie, Gabriel B. Kwarko, Emmanuel Broni, Michael B. Adinortey and Michael D. Wilson

Submitted: August 10th, 2021 Reviewed: September 23rd, 2021 Published: December 18th, 2021

DOI: 10.5772/intechopen.100594

Chapter metrics overview

133 Chapter Downloads

View Full Metrics


Trypanothione reductase (TR), a flavoprotein oxidoreductase is an important therapeutic target for leishmaniasis. Ligand-based pharmacophore modelling and molecular docking were used to predict selective inhibitors against TR. Homology modelling was employed to generate a three-dimensional structure of Leishmania major trypanothione reductase (LmTR). A pharmacophore model used to screen a natural compound library generated 42 hits, which were docked against the LmTR protein. Compounds with lower binding energies were evaluated via in silico pharmacological profiling and bioactivity. Four compounds emerged as potential leads comprising Karatavicinol (7-[(2E,6E,10S)-10,11-dihydroxy-3,7,11-trimethyldodeca-2,6-dienoxy]chromen-2-one), Marmin (7-[(E,6R)-6,7-dihydroxy-3,7-dimethyloct-2-enoxy]chromen-2-one), Colladonin (7-[[(4aS)-6-hydroxy-5,5,8a-trimethyl-2-methylidene-3,4,4a,6,7,8-hexahydro-1H-naphthalen-1-yl]methoxy]chromen-2-one), and Pectachol (7-[(6-hydroxy-5,5,8a-trimethyl-2-methylidene-3,4,4a,6,7,8-hexahydro-1H-naphthalen-1-yl)methoxy]-6,8-dimethoxychromen-2-one) with good binding energies of −9.4, −9.3, 8.8, and −8.5 kcal/mol, respectively. These compounds bound effectively to the FAD domain of the protein with some critical residues including Asp35, Thr51, Lys61, Tyr198, and Asp327. Furthermore, molecular dynamics simulations and molecular mechanics Poisson-Boltzmann surface area (MMPBSA) computations corroborated their strong binding. The compounds were also predicted to possess anti-leishmanial activity. The molecules serves as templates for the design of potential drug candidates and can be evaluated in vitro with optimistic results in producing plausible attenuating infectivity in macrophages.


  • Leishmania
  • trypanothione reductase
  • oxidative stress
  • natural product
  • pharmacophore modeling
  • virtual screening
  • molecular dynamics

1. Introduction

Leishmaniasis is a disease caused by a single-cell eukaryotic parasite of the Leishmaniaspecies. This protozoan parasite causes a substantial level of morbidity and mortality. Leishmaniahas a digenetic life cycle [1]. In mammals, the parasite colonizes macrophages, transforming into intracellular amastigotes. The parasite has an adaptive way to life conditions. The amastigotes tolerate low pH and are hydrolase resistant [2].

Trypanothione is a major product of the trypanothione biosynthesis pathway in trypanosomes which is crucial in maintaining cellular redox potential and is essential for the parasite’s survival. This molecule is catalyzed by so many enzymes for which Leishmania majortrypanothione reductase (LmTR, E.C. plays a critical role in the biosynthetic pathway. TR reduces trypanothione (T[S]2) to dithiol (T[SH]2). They catalyze the transfer of electrons from NADPH to their specific substrate via an FAD prosthetic group [3]. The reduced form is critical in regulating oxidative stress by reacting with reactive oxygen species (ROS) that are produced by the macrophage. T[SH]2 is not only needed for detoxification of peroxides but also required for the synthesis of DNA precursors, homeostasis of ascorbate, sequestration and export of thiol conjugate [4].

Trypanothione reductase is a member of the disulphide oxidoreductase family of enzymes. It has an analogue in the human body, glutathione reductase (GR) which also carries out oxidoreductive reactions. But LmTR does not process GSSG and host GR does not reduce T[S]2 [5, 6]. The ascribed reasons for targeting LmTR include the following: (i) trypanothione reductase is a key enzyme in regulating a reducing environment aiding in disease pathogenesis, (ii) this parasite does not depend on the host for reduced trypanothione, (iii) it has a less close known homologous protein in humans; (iv) the availability of template homologs for modeling purpose; and (v) moreover, in vitrotrials have proven LmTR to be a good therapeutic target [7].

Several inhibitors have been screened against this enzyme causing a reduction of infectivity and decreased capacity of the parasite to survive within intracellular macrophages. Potent compounds, such as 7-chloro-4-nitro-5-quinazolin-4-ylsulfanyl-2,1,3-benzothiadia-zole (CNQB) and 4-phenyl-5-(4-nitro-cinnamoyl)-1,3,4-thiadiazolium-2-phenylamine-chloride (PNTPC) with IC50 values 0.58 and 1.63 μM, respectively have already been tested in an in vitroassay against trypanothione reductase of trypanosomatids [8, 9].

Computer-aided drug designing is an in-silicoapproach for drug discovery that combines computational and pharmaceutical research [10]. This application helps in spanning the drug discovery pipeline and helps to speed up and rationalize the drug design process while reducing costs [11]. Ehrlich in 1909 first defined the term pharmacophore as ‘a molecular framework that carries (phoros) the essential features responsible for a drug’s (pharmacon) biological activity [12]. These features are essential functional groups of atoms in a three-dimensional position that interact with a receptor. Ligand-based drug design can be performed in association with molecular docking. These methods can be combined to identify a number of new hit compounds with potent inhibitory activity and to understand the main interactions at the binding sites. Appropriate use of these methods can improve the ability to identify and optimize hits and confirm their potential to serve as scaffolds for producing new therapeutic agents [13].

Drugs currently used for the treatment of human leishmaniasis are toxic, having severe adverse reactions which limit their use. Aside this includes, increase in resistance by the parasite, high cost of available drugs, lack of efficacy against VL\HIV co-infections with standard chemotherapy, and the development of a single drug or formulation for all forms of leishmaniasis [14, 15, 16, 17]. Therefore, the development of novel, effective drugs with reduced side effects, is still a major priority for health researchers, in spite of many compelling research reports published on antileishmanial agents in the last 10 years [18]. In this study, in silicomethod of identifying leads was used incorporating the knowledge of pharmacophore and virtual screening to arrive at lead molecules. The overall goal of the study was to predict with high degree potent selective inhibitors of LmTR from the African Natural Product Database (AfroDb) and North African Natural Product Database (NANPDB).


2. Materials and methods

2.1 Protein homology modeling

The protein sequence of trypanothione reductase of Leishmania majorwas retrieved from the NCBI database with the accession number XP_001687512.1, having 491 amino acids [19]. Using the Basic Local Alignment Search Tool (BLAST), the query sequence was compared to known structures which generated structures similar to LmTR. Modeller 9.20 [20] was used for modeling the structure of LmTR.

2.2 Active site prediction and quality assessment

To predict the active site, the protein was submitted to CASTp 3.0 [21, 22]. The predicted active site was corroborated via a blind docking process using AutoDock Vina within PyRx version 0.9.7 [23, 24]. The active region was confirmed by the ‘Toggle selection of Spheres’ function which highlighted predicted residues from CASTp 3.0. Binding pocket was also viewed with PyMOL v2.0.0 [25]. The quality of the modeled protein was assessed by some quality measure tools. This included PROSA which determines the quality of experimentally solved structures and theoretical models in protein engineering by comparing these to that of experimentally solved protein structures in the PDB database [26]. Verify3D was used to validate the three-dimensional structure of the model [27]. PROCHECK, a quality assessment tool was also used to check the stereochemical properties of the model by generating a Ramachandran plot [28]. ProQ was also used to carry out further validation. ProQ predicted protein quality based on the LGscore and MaxSub scores [29].

2.3 Energy minimization of protein target

The modeled LmTR was energy minimized using GROMACS 5.1.1 [30, 31]. Simulations were performed with the force field, GROMOS96 43a1. The system was solvated using an equilibrated SPCE216 water model. The charged protein had a net charge of −9 which was neutralized with an equal amount of Na+ ions. Energy minimization of the protein was then carried out.

2.4 Pharmacophore modeling and screening

2.4.1 Pharmacophore generation

Pharmacophore model, virtual screening, and molecular docking studies were performed to find novel LmTR inhibitors. The ligand-based structural design incorporates the absence of the macromolecular structure by generating pharmacophore models from a set of ligands. This method takes advantage of the conformational flexibility of the ligands [32]. The active compounds, CNQB and PNTPC were retrieved with IDs CID1323435 and CHEMBL242165 from PubChem ( CHEMBL ( [33, 34], respectively and were used to generate the pharmacophore. These ligands served as a training set for pharmacophore generation. All customized settings were kept in default.

2.4.2 Library preparation and pharmacophore validation

LigandScout 4.3 was used for screening a total of 5813 compounds including 885 AfroDb entries found in ZINC database [35] and 4928 NANPDB compounds [36]. The two actives were converted to SMILES format and submitted to the Directory of useful decoys and enhanced (DUD-E) database to generate decoys for the screening [37]. A total of 100 decoys were generated and used as a decoy library. The libraries were then converted into a .ldb file format. The reliability of the pharmacophore model was validated by the area under the receiver operating characteristic (ROC) curve (AUC) [38] using two descriptors, selectivity and sensitivity.

2.4.3 Screening with pharmacophore

LigandScout 4.3 allows in-silicoscreening of compound libraries using pharmacophore models as filter criteria. Database of active ligands was selected and marked in green with decoys marked in red. The screening process was initiated to generate hits corresponding to the pharmacophore model. These compounds were saved in an “sdf” file format to be used in the docking process.

2.5 Molecular docking

2.5.1 Docking of LmTR

The generated hit compounds were uploaded into PyRx (Version 0.9.6) [23]. The energy of the ligands was minimized using Universal Force Field (UFF) option in Open Babel incorporated in PyRx prior to docking. This was done to obtain 3D ligand structures which constitute atomic elements that have proper bond lengths between their atoms [39]. Ligands were converted to PDBQT format using AutoDock Vina embedded in the PyRx. Predicted active site residues were selected within a grid box of dimensions X: 42.39 Å, Y: 35.47 Å, and Z: 31.05.14 Å; and centre X: 28.58 Å, Y: 57.09 Å, and Z: −2.24 Å within the AutoDock Vina environment of PyRx for docking process.

2.5.2 Docking validation with AUC

Validation of the algorithm used for the docking process was carried out by generating an AUC plot. Decoys of five known inhibitors in complex with Leishmaniatrypanothione reductase of several species of Leishmaniawhich included 2,3,4,6-tetra-O-acetyl-1-thio-beta-beta-d-glucopyranose (Auranofin); 4-[[1-(4-ethylphenyl)-2-methyl-5-(4-methylsulfonylphenyl)pyrrol-3-yl]methyl]thiomorpholine (CHEMBL1277380); 6-sec-butoxy-2-[(3-chlorophenyl)sulfanyl]-4-pyramidinamine (RDS); 2-(diethylamino)ethyl4-((3-(4-nitrophenyl)-3-oxopropyl)amino)benzoate (ZINC8782981); {N}-(4-azanylbutyl)- ∼ {N}-(2-azanyl-2-oxidanylidene-ethyl)-7-(3-azanyl-3-oxidanylidene-propyl)-4-(dimethylamino)-2-(2-naphthal en-2-ylethylamino)pyrrolo[2,3-d]pyrimidine-6-carboxamide (H6H) were retrieved from DUD-E to generate an AUC curve. The result was correlated and plotted using their respective binding energies as the only variable via easyROC (Ver. 1.3) [40].

2.5.3 Docking validation via superimposition and alignment

Superimposition of the crystallographic ligand and re-docking poses was used as a means of validating docking. The five crystallographic ligands in complex with Leishmaniatrypanothione reductase were removed from the co-crystallized complex and re-docked. The pdb files of the re-docked complex were uploaded into PyMOL together with the solved complex from the Protein Data Bank. LigAlign [41] was then used to calculate the root mean square deviation between the superimposed re-docked and co-crystalized ligands. Superimposition also allowed the identification of critical overlapping residues via LigPlot+. The FAD molecule from the template selected for modeling (PDB ID: 2JK6) was also extracted and docked in LmTR. A comparative study of the original FAD-2JK6 complex and FAD-LmTR complex was done using LigPlot+.

2.6 Identification of lead compounds

To identify lead compounds the binding energy, molecular bond interactions, pharmacological, and physiochemical properties were considered. This step helped to filter generalized hit compounds. These compounds in SMILES format were submitted to SwissADME [42], which calculates the corresponding ADME (absorption, distribution, metabolism, and excretion) properties of the compounds. Hydrogen bond interactions of the ligand–protein complexes were studied using LigPlot+ and PyMOL.

The hit compounds were physiochemically profiled to identify their drug-likeness and solubility in water. Lipinski’s rule of 5 was used as a metric to narrow down druggable compounds [43]. Pharmacokinetic properties of predicted compounds were determined in silico. This included cytochrome inhibition, P-glycoprotein (P-gp) substrates, gastrointestinal (GI) absorption, and the blood–brain barrier (BBB) permeant.

2.7 Prediction of activity spectra for substances (PASS) for leads

PASS assesses the probability that a compound has a suspected biological activity [44]. It has been well known that each substance has a wide spectrum of biological activities as evident from some new uses of many old drugs. The SMILES format of the leads were submitted to [45] to predict possible biological activity.

2.8 Molecular dynamics simulation and MM-PBSA calculation of protein–ligand complexes

By employing GROMACS 2018 [31], the chain A of LmTR and the LmTR–ligand complexes were subjected to molecular dynamics simulations. Ligand topologies were generated via PRODRG which were converted to .gro files. Solvation of each of the protein–ligand complexes in a dodecahedron box was followed by neutralization of the output with sodium and chlorine ions. Each complex was minimized using the steepest descent algorithm coupled with the GROMOS43A force field. Equilibration protocol was carried out to restrain and relax protein–ligand positions. The MD production run was carried out for 100 ns. The output file was used in downstream processes to generate root mean-square deviation (RMSD), root mean-square fluctuation (RMSF) and radius of gyration (Rg) plots with Xmgrace ( Molecular mechanics Poisson–Boltzmann surface area (MM-PBSA) was employed in calculating the free energies of the complexes. MM-PBSA was carried out using g_mmpbsa, which calculates binding energy components and the per-residue energy decomposition [46]. Graphs were generated using R-programming showing energy interactions.


3. Results and discussion

3.1 Homology modeling

L. majortrypanothione reductase was modeled with an appreciable degree of accuracy and validated via several bioinformatics tools (Figure 1A). The protein sequence of L. majorwith the accession number XP_001687512.1, having 491 amino acids when searched amongst protein homologs in NCBI resulted in a list of similar structures with PDB codes 2JK6, 2X50, 6ER5, 1FEA from Leishmania infantumwith their respective percentage identity of 95.72, 95.71, 95.70, and 78.78%, respectively. The similarity of structure and sequence between LmTR protein sequence and the template 2JK6 was 95.72%. This favored its selection for the modeling process. The best model selection was based on the least discrete optimization potential energy (DOPE) score which informs about the energy of the protein.

Figure 1.

(a) 3D structure of modeled protein in cartoon representation. A monomer ofLmTR with alpha helices represented in red, the beta sheets in yellow and the loops in green. (b) Surface representation of theLmTR. Active site (in red) is the FAD binding region.

3.2 Active site prediction

In predicting the active site of the protein, results obtained showed 73 binding pockets. Several pockets were identified but the pocket with the largest volume and surface area of 595.278 Å3 and 924.887 Å2, respectively was selected. Larger pockets favor conformational rotation during virtual screening. A total of 80 amino acid active site residues were predicted from CASTp 3.0 (Table A1). The predicted binding site was visualized using PyMOL (Figure 1B).

The predicted active site was finally confirmed to be FAD binding site. TR has been studied to be a homodimer protein constituting FAD-binding, NADPH-binding, central, and interface domains [5]. The predicted site was concluded to be an FAD site in LmTR by comparative studies between the FAD-2JK6 and FAD-LmTR complexes (Figure A1-A). The study showed common hydrogen bonding residues including Ser14, Gly15, Arg287, and Thr335 interacting with FAD. These residues correlate with that which was predicted by CASTp 3.0 including other several common residues participating in hydrophobic interactions such as Gly13, Asp35, Ala46, Gly50, Thr51, and Cys52 (Figure A1-A).

3.3 Structural validation and quality prediction

PROSA was used to determine the quality of the structure by comparing it to protein structures that are experimentally solved in the PDB database. It validated the model based on the “quality score or z-score” with a value of −11.68. The z-score shows whether the predicted model is of X-ray or NMR quality with regards to the amino acid residues length. This z-score value showed that the modeled protein fell in the range of proteins solved experimentally by X-ray crystallography (Figure A2-A). Verify3D validated the 3D structure of the model. A good 3D structure is expected to have at least 80% of its amino acids to have scored greater than or equal to 0.2 in a 3D/1D profile. This model passed with an appreciable result of 91.65% of residues having a score greater or equal to 0.2 (Figure A2-B). Further validation with PROCHECK resulted in the generation of a Ramachandran plot (Figure A3). The plot described the rotations of the polypeptide backbone around the bonds between N-Cα (Phi, φ) and Cα-C (Psi, ψ). The plot allowed the viewing of the distribution of these torsional angles taking into consideration the allowed rotations and rotations that are unfavored which can result in a collision or steric hindrance. A protein with 90% of its residues in the most favorable region is considered a good model. ProQ predicted the quality of protein based on the MaxSub and LGscore scores [29]. The predicted LG score was 6.447 and MaxSub score was 0.520. LGscore >4 implied that the model was extremely good. Also, a MaxSub score > 0.5 implied that the model was very good [47].

3.4 Pharmacophore modeling

The active ligands used as training sets to develop a pharmacophore allowed features similar to the two compounds to be identified and combined into a single geometric function as the basis for the generation of the pharmacophore (Figure 2). Pharmacophore generated utilized features that contributed regions of hydrophobicity and hydrogen bond acceptors incorporated in the model for selective screening. The oxygen from nitrogen dioxide contributed to hydrogen bond acceptors with aromatic and alkene groups contributing to the hydrophobic region (Figure A4). A number of 10 hypotheses were developed for the model. The best hypothesis with a similarity of 58.12% had a score of 0.8537 and was selected based on the AUC score of 0.99 generated by LigandScout (Figure A5-A). The AUC score was used as a metric to validate how best the pharmacophore model created could distinguish rightly between active compounds and decoys. This intends to reduce false positives and negatives during the screening process.

Figure 2.

2D structure of PNTPC (left) and CNQB (right) showing shared features. PNTPC and CNQB showed two hydrogen bond acceptors. Both compounds had an aromatic ring as a common feature. These served as active ligand for the training set.

The screening process was successfully completed with the model which identified 42 compounds that matched the pharmacophore model with a pharmacophore fit score ranging from 55.32 to 57.98 (Table A2).

3.5 Docking validation

The docking system was validated using the five inhibitors and 250 decoys generated with DUD-E and further used their binding energies to plot an AUC with a value of 0.702 (Figure A5-B). AUC value of 1.0 verifies that the prediction of hits obtained from the hypothesis is perfect whereas values of 0.5 and less than 0.70 imply average and moderate random selection respectively [48].

Furthermore, the validation of the molecular docking was also undertaken by aligning the re-docked ligands with their respective co-crystallized complexes taken from PDB. The RMSD values of the alignment between the re-docked ligands and the co-crystalized ligands in complexes of 2YAU, 4APN, 5EBK, 6ER5, 6I7N were 1.483, 3.020, 1.920, 2.712, and 2.465 Å, respectively. Only two of the RMSD values (5EBK and 2YAU) of the alignments were below 2.0 Å, which is considered the threshold for good alignment.

Superimposition also validated the accuracy of docking at the predicted active site. The FAD molecule from the template selected for modeling, 2JK6 when extracted and docked in LmTR, showed similar surrounding residues in the pocket of FAD-2JK6 and FAD-LmTR complex. Ligand alignment of these two complexes gave an RMSD of 3.291 Å which is above the expected threshold. Nonetheless, the FAD-LmTR docking simulated a pose that showed common residues such as Ser14, Gly15, Arg287, and Thr335 taking part in hydrogen bonding in these complexes (Figure A1). All these verified that the docking system performed very well in docking ligands to the active site.

3.6 Virtual screening of pharmacophore hits

When docking validation was verified, molecular docking was carried out. Molecular docking predicted various conformations of each ligand in the binding site of the LmTR. Compounds were selected based on their binding energies. The binding energies gave a theoretical value that relates the affinity of the ligand to the protein model. The results of the respective pharmacophore fit scores and binding energies are well documented (Table A2). The 42 compounds obtained were narrowed down to 11 by considering their binding energies and pharmacophore fit scores (Table 1). With respect to these studies to find inhibitors of the FAD binding site for which FAD is a prosthetic group that already binds tightly to the catalytic site, it was reasonable to select the compounds with binding energies below −9.0 kcal/mol or relatively closer to −9.0 kcal/mol which can provide a plausible binding advantage when acting as inhibitors at the site. The compound ZINC95486081 had the lowest binding energy of −9.8 kcal/mol and the highest was observed from evoxine as −6.2 kcal/mol (Table A2).

Predicted ligandsPharmacophore fit scoreBinding energy/(kcal/mol)Hydrogen bond Residues and length (Å)Hydrophobic bond interacting residues
ZINC9548608155.95−9.8Lys60 (2.92)
Ser178 (2.87)
Gly13, Gly15, Asp32, Ala46, Thr51, Cys52, Val55, Gly56, Ala159, Thr160, Tyr198, Arg287, Asp327, Met333, Leu334, Thr335
MTPA56.37−9.4Thr51 (3.12)
Thr293 (2.89)
Asp327 (3.13)
Ser14 (3.17)
Gly11, Gly13, Gly15, Asp32, Val36, Ala46, Gly50, Cys52, Val55, Gly56, Ala159, Thr160, Tyr198, Arg287, Met333, Leu334, Thr335
Karatavicinol56.5−9.4Thr51 (3.28)
Arg287 (3.00)
Thr293 (3.01, 2.97)
Asp327 (2.82)
Ser14, Gly11, Gly13, Val34, Asp35, Val36, Gly50, Cys52, Ala46, Phe126, Gly127, Ala159, Thr160, Gly161, Tyr198, Arg290, Met333, Ala338
Taccalin56.42−9.4Ser14 (3.21)
Ala365 (3.25)
Gly13, Thr51, Cys52 Gly56, Cys57, Lys60, Gly161, Ile199, Thr198, Gly326, Met333, Leu334, Thr335, Ala338
Marmin56.18−9.3Val34 (3.17)
Thr51 (2.99, 3.04)
Thr160 (2.83)
Thr335 (2.87)
Leu10, Gly11, Gly13, Ser14, Asp35, Ala46 Gly50, Cys52, Gly127, Ala159, Gly161, Arg290, Leu294, Ala327, Leu334, Ala338
13-Hydroxyfeselol55.62−9.1Val362 (2.79)
Thr374 (2.85)
Gly376 (3.20, 3.26)
Lys60, Thr198, Gly229, Phe230, Gly326, Leu334, Cys364, Ala365
Betaxanthin57.51−8.9Ser14 (2.85)
Cys52 (3.03)
Gly127 (3.03, 3.22, 3.30)
Gly11, Gly13, Val34, Asp35, Val36, Gly50, Cys52, Ala46, Phe126, Ala159, Thr160, Gly161, Tyr198, Arg290, Asp327 Met333, Ala338
Colladonin55.90−8.8Asn330 (2.85)Lys60, Gly197, Tyr198, Tyr221, Arg287, Phe230, Leu334, Ala365
Feselol55.63−8.8Asn330 (2.87)Lys60, Gly197, Tyr198, Tyr221, Arg287, Phe230, Leu334, Ala365
ZINC3865803555.95−8.7Tyr198 (3.31)
Val362 (2.92)
Tyr374 (3.30)
Gly376 (2.86, 3.05)
Ile199, Phe230, Gly286 Arg287, Met333, Leu334, Cys364, Cys375
Pectachol57.18−8.5Lys60 (2.93)
Gly376 (3.11)
Tyr198, Gly229, Phe230, Val332, Met333, Leu334, Ala365, Val362, Cys364, Val366, Phe367
FAD molecule and inhibitors from 6ER5 and 4APN
ZINC8782981_−7.2Lys60 (3.12)
Arg287 (3.08)
Cys52, Cys57, Tyr198, Gly229, Phe230, Val332, Met333, Leu334, Ala365, Val362, Cys364
CHEMBL1277380_−8.2Lys60 (2.89)Tyr198, Gly229, Phe230, Val332, Met333, Leu334, Ala365, Val362, Gly376
FAD_−9.0Ser14 (3.11)
Gly15 (2.96)
Ala159 (2.97)
Tyr198 (2.87)
Arg287 (2.75)
Met333 (2.78)
Thr335 (2.88, 3.31)
Gly13, Gly50, Thr51, Cys52, Ser162, Gly197, Gly229, Phe230, Asp327, Leu334, Ala338

Table 1.

The table shows parameters involved in the selection of lead compounds. This included pharmacophore fit score, the binding energy and the number of hydrogen and hydrophobic bond interacting residues.

3.7 Protein–ligand interaction

Molecular interaction studies are important for understanding the mechanism of biological regulation at the molecular level and as such also provides a theoretical basis for drug design and discovery [49, 50]. Hydrogen and hydrophobic interactions are key players in stabilizing energetically favored ligands, in an open conformational environment of protein structures [29]. The intermolecular interaction and bond lengths of these 11 compounds were observed. The compound which showed the highest binding affinity, ZINC95486081 formed two hydrogen bonds with residues Lys60 and Ser178 with respective bond lengths of 2.92 Å and 2.87 Å. Five compounds including MTPA, Karatavicinol, Marmin, Betaxanthin and ZINC38658035 had four residues as the highest number of residues partaking in hydrogen bonding. MTPA formed hydrogen bonds with residues Thr51, Thr293, Asp327, and Ser14. Karatavicinol on the other hand formed hydrogen bonds with Thr51, Arg287, Thr293, and Asp327 (Figure A6-A). Marmin formed hydrogen bonds with Val34, Thr51, Thr160, and Thr335 (Figure A6-B). Betaxanthin also bonded with Ser14, Cys52, Gly127, and Thr335. Finally, ZINC38658035 formed hydrogen bonds with Tyr198, Val362, Tyr374, and Gly376. The compound 13-hydroxyfeselol was the only hit that formed three hydrogen bonds with Val362, Thr374, and Gly376. Taccalin and Pectachol formed only two hydrogen bonds with Ser14 and Ala365. On the other hand, Colladonin and Feselol formed the least hydrogen bond residues with Asn330. The shortest bond length of 2.83 Å was exhibited by Marmin with Thr160. Betaxanthin showed the highest pharmacophore fit score of 57.51 followed by Pectachol (57.18). 13-Hydroxyfeselol showed the lowest fit score of 55.62.

3.8 Pharmacological profiling

To identify lead compounds, the binding energy, molecular bond interactions, pharmacological, and physiochemical properties were considered. This step helped to filter generalized hit compounds. The top 11 compounds were profiled in silicoto characterize compounds with drug-likeness and good water solubility. Lipinski’s rule of 5 was used as a metric to narrow down druggable compounds. The rule factors in the compound’s molecular weight which should not exceed 500 g/mol, hydrogen bond donors must not be more than 5, log-p value must be less than or equal to 5 and hydrogen bond acceptors must not be more than 10. All the 11 top hits passed Lipinski’s rule of 5 with a good bioavailability score of 0.55 (Table A3).

3.9 Pharmacokinetic properties

Further filtering analysis subjected all 11 pharmacophore hits to pharmacokinetics profiling taking into consideration parameters such as gastrointestinal (GI) absorption, blood–brain barrier (BBB) permeation, the permeability of glycoprotein (Pgp), and cytochrome P450 (CYP). Physical parameters such as drug solubility may affect oral bioavailability but in most cases, the major determining factors are likely to be metabolism by CYP and absorption at the intestinal level [51]. CYP3A4 has been known to be responsible for the metabolism of about 50% of all drugs [52] and therefore inhibition of cytochrome can affect oxidation of substrates in cells. Absorption of drugs in the intestine if found high favors the efficacy of the compound as a drug. Multi-drug resistance transporters, such as P-glycoproteins, are essential for many cellular processes that require the transport of substrates across cell membranes [53]. Compounds that are P-gp substrates may face continual efflux which can affect the efficacy of drugs. The blood–brain barrier (BBB) prevents the brain uptake of most pharmaceuticals [54]. This is a disadvantage to neurological diseases but would be of merit since the disease of study is not related to the brain. Compounds that cross the blood–brain barrier may elucidate unwanted biological activities that could be dangerous to health. Therefore, the negative inference would be good for the compound. ZINC95486081 was predicted to show inhibition to three CYP isoenzymes. Karatavicinol, ZINC38658035, and Marmin excelled with an appreciable result (Table A4). For the purpose of narrowing down leads with potential for further computational analysis, compounds with low gastrointestinal absorption were side-lined. This included Taccalin and Betaxanthin.

3.10 Prediction activity spectra for substance (PASS)

The biological activity of the selected drug-like candidates was then evaluated using PASS. It is well known that each substance has a wide spectrum of biological activities as evident from some new uses of many old drugs. This allows the tool to utilize this information to predict biological activities based on their probable activity (Pa) and probable inactivity (Pi). When Pa is greater than Pi (Pa > Pi), the compound is likely to possess the predicted biological activity [55, 56]. PASS predicted Karatavicinol, Marmin, Colladonin and Pectachol to be potential antileishmanial agents (Table A5). Colladonin showed the highest Pa of 0.768 and Pi of 0.006 followed by Taccalin (Pa of 0.711 and Pi of 0.009) and Pectachol with a Pa of 0.694 and Pi of 0.009. Betaxanthin had no prediction as an antileishmanial agent.

3.11 Selection of lead compounds for MD and MM-PBSA analysis

The various lead compounds were considered for selection based on the criteria above. ZINC95486081 and MTPA compound although had high binding energy trailed in pharmacokinetic properties and showed Pa less than 0.500. We eliminated Taccalin and Betaxanthin because of their low GI absorption and low Pa values (Tables A3 and A4). Compounds predicted with good probable activity for antileishmanial activities included Karatavicinol, Taccalin, Marmin, 13-hydroxyfeselol, Colladonin, Feselol and Pectachol. A literature search revealed Feselol to have antiprotozoal activity against Trypanosoma brucei(IC50 8.1 μM), Trypanosoma cruzi(IC50 8.6 μM), and Leishmania infantum(IC50 6.8 μM). Similarly, 10′R-karatavacinol has presented activity against T. brucei(IC50 32.4 μM), T. cruzi(IC50 9.4 μM), and L. infantum(IC50 32.4 μM) [57]. That of Feselol and hydroxyfeselol was eliminated because it had been worked on experimentally and also extracts from Ferula genuswhich is known to exhibit antiviral, antibacterial, and antileishmanial properties [58]. Marmin and Pectachol presented optimistic potential to look into. The high number of hydrogen bonds and binding energy allowed for the inclusion of Karatavicinol in the MM-PBSA to observe their stabilization with this protein target. Also, the high prediction of Colladonin as a compound with probable activity also increased the chance of this compound for MM-PBSA analysis. Therefore, Marmin Pectachol, Karatavicinol, and Colladonin were considered for further analysis in molecular dynamics.

3.12 Molecular dynamic simulation of protein–ligand complex

Molecular dynamics simulation allowed the early view of proteins as relatively rigid structures to be replaced by a dynamic model in which the internal motions and resulting conformational changes play an essential role in its function [59]. An RMSD plot generated after molecular dynamics simulation showed a deviation of about 0.25 Å (Figure A7). Further scrutinized with molecular dynamics simulations gave the protein a dynamic dimension to its 3D structural form producing a realistic environment for the ligand interactions that were carried out in the docking process.

Molecular dynamics simulations can also capture a wide variety of important biomolecular processes, including conformational change, ligand binding, and protein folding [60]. The stability of docked protein–ligand complexes was determined by their (RMSD) plots generated from the MD simulation output file. The backbones of the four complexes were observed to be stable over time (Figure 3). The fluctuations of the protein–ligand complexes were analyzed within the system to check for movement and structural stability during the course of the simulation. These movements and stability are significant for the complex functioning inside living systems. The backbone of the LmTR-ZINC8782981 complex showed the greatest stability with an average RMSD of 0.25 nm amongst all the complexes. The LmTR-CHEMBL1277380 complex was fairly stable with RMSD of 0.4 nm. Amongst the four leads, Pectachol was found with the lowest RMSD of 0.3 nm. In terms of stability, the compound Marmin and Colladonin proved to be very much stable around 0.73 nm over the production time of 100 ns. The RMSD of Karatavicinol was observed to show stability from 0 to 60 ns and increased its RMSD to 0.5 nm from 60 to 100 ns.

Figure 3.

RMSD values of theLmTR-ligand complexes of the four leads (Karatavicinol, Marmin, Pectachol, and Colladonin) and the two known inhibitors after 100 ns. The complexes in the graph are color coded.

The flexibility of residues contribution by the LmTR was assessed by the root mean square fluctuation (RMSF). RMSF indicates the flexibility of different regions of a protein, which can be related to crystallographic B factors [61]. The results of the RMSF plots showed consistency for the docked complexes (Figure A8). The highest fluctuations exhibited was observed around residue numbers 70–90, with Karatavicinol and CHEMBL1277380 showing higher fluctuation levels followed by ZINC8782981 and Marmin complexes. Other regions where good fluctuations were observed include residues between 395 and 410 and 465–480. Overall, Marmin showed more fluctuations around most residues with a distinct difference at residue numbers 260 and 310.

The compactness of the complexes over simulation time is determined by the Rg. If a protein is folded well, it will likely maintain a relatively steady value of Rg, whereas its value will change over time if the protein unfolds [62]. Rg values of all complexes indicated stable complexes over 100 ns (Figure A9). The Rg graph showed most compounds experienced a fairly stable Rg. Marmin experienced the lowest Rg value around 2.33 nm compared to other complexes. This was followed by Colladonin, Pectachol and Karatavicinol with Rg values of 2.37, 2.42, and 2.45 nm, respectively. Between the known inhibitors, CHEMBL1277380 was observed to have an average Rg value of 2.46 nm whilst ZINC8782981 showed the average highest value of around 2.5 nm. Inferring from the Rg graph, the compactness of the LmTR–Marmin, –Colladonin and –Pectachol complexes were maintained after complex formation.

3.13 Evaluation of leads using MM-PBSA

MM-PBSA was employed to calculate free binding energies by per-residue decomposition of the protein complexes. At a quantitative level, simulation-based methods provide substantially more accurate estimates of ligand binding affinities (free energies) than other computational approaches such as docking [63]. Residues contributing binding free energy greater than 5 kJ/mol or less than −5 kJ/mol are considered critical for binding of a ligand to a protein [64]. MM-PBSA results showed only Asp327 amongst the hydrogen bonding residues of Karatavicinol to contribute a per residue decomposition energy of 13.65 kJ/mol. Amino acid residue Asp35 (21.89 kJ/mol) was observed with such greater contribution (Figure A10). The complex of LmTR–Marmin also showed surrounding hydrophobic residues Asp35 (−8.62 kJ/mol), Ala46 (−7.65 kJ/mol), Arg290 (8.83 kJ/mol), and Glu141 (−5.12 kJ/mol) with their energy decomposition to be greater than 5 kJ/mol and less than −5 kJ/mol. The only hydrogen bonding residue that showed a relevant contribution of energy decomposition was Thr51 (−16.36 kJ/mol) (Figure 4). Hydrophobic amino acid residues Lys61 (−5.16 kJ/mol), Tyr198 (−11.29 kJ/mol), Asp327 (5.56 kJ/mol), and Arg331 (6.38 kJ/mol) showed relevant contribution to the total binding energy of the LmTR–Colladonin complex (Figure A11). Moreover, only hydrogen bonding residue Lys60 (13.32 kJ/mol) in LmTR–Pectachol complex showed to be a critical residue in binding. Other surrounding residues contributed substantially to the per residue energy decomposition in the LmTR–Pectachol complex. This included Lys61 (−10.97 kJ/mol), Arg287 (−6.91 kJ/mol), Asp327 (9.30 kJ/mol), Met333 (−6.72 kJ/mol), Leu334 (−8.91 kJ/mol), Lys361 (6.46 kJ/mol), and Cys364 (−6.07 kJ/mol) (Figure A12). Deducing from the substantial contribution of energy per decomposition of residues, we propose Asp35, Thr51, Lys61, Tyr198, and Asp327 to be critical in intermolecular bonding and stabilization of ligands at the FAD active site.

Figure 4.

MM-PBSA plot of the binding free energy decomposition contribution per residue ofLmTR–Marmin complex. Coded red lines represent surrounding active site amino acid residues.

3.14 Other energy terms

Van der Waals forces, electrostatic and polar solvation energies, and SASA are relevant energy terms contributing to the overall free binding energy of the complex. The van der Waals energy refers to the weak attraction existing between the intermolecular forces. The van der Waals energy observed in our study showed Karatavicinol and CHEMBL1277380 to have the lowest and highest energy of −228.565 and − 171.823 kJ/mol, respectively. Colladonin, Marmin, and Pectachol also showed relatively low van der Waals energy of −189.289, −189.229, and − 209.538 kJ/mol, respectively as compared with ZINC8782981 with −222.123 kJ/mol. Electrostatic energy refers to the potential energy of a system consisting of different electric charges [65, 66, 67]. The lowest electrostatic energy was exhibited by Marmin (−386.401 kJ/mol) followed by Pectachol (−286.260 kJ/mol), and Colladonin (−249.067 kJ/mol). Karatavicinol and the other two inhibitors were observed with high electrostatic energy (Table 2). Some studies have observed that van der Waals and electrostatic forces contribute favorably to the energetics of binding along with simulations that favor the binding of complexes [66, 68].

CompoundVan der Waals energy (kJ/mol)Electrostatic energy (kJ/mol)Polar solvation energy (kJ/mol)SASA energy (kJ/mol)Binding energy (kJ/mol)
ZINC8782981−222.123 ± 14.568−52.495 ± 17.662245.049 ± 33.654−24.829 ± 0.962−54.399 ± 20.084
CHEMBL1277380−171.823 ± 14.173−2.926 ± 5.48582.884 ± 16.706−19.869 ± 1.202−111.732 ± 16.514
Karatavicinol−228.565 ± 12.673−32.345 ± 21.415227.483 ± 27.305−24.217 ± 1.100−57.644 ± 24.019
Marmin−189.289 ± 16.726−386.401 ± 30.540484.074 ± 28.991−17.498 ± 1.050−109.114 ± 23.461
Pectachol−189.229 ± 18.203−286.260 ± 49.152430.604 ± 71.136−18.602 ± 1.308−63.487 ± 33.289`
Colladonin−209.538 ± 18.908−249.067 ± 40.851427.216 ± 49.348−17.548 ± 1.122−48.936 ± 24.773

Table 2.

The energy terms obtained after MM-PBSA analysis of the protein–ligand complexes.

The energy values are presented as mean ± standard deviation kJ/mol.

Polar solvation energy on the other hand refers to the electrostatic interaction that exists between the solute and the continuum solvent [69]. The highest polar solvation energy amongst the leads was exhibited by Marmin (484.074 kJ/mol) and the lowest by Karatavicinol (227.483 kJ/mol). Solvent accessible surface area (SASA) energy was calculated after MD. This represents the non-polar solvation energy [69]. This energy measures the interactions that exist between the complex and the solvents. Amongst the leads, Karatavicinol obtained the lowest SASA energy followed by Pectachol, Colladonin, and Marmin (Table 2). Relative to these were the low SASA energies of the inhibitors ZINC8782981 and CHEMBL1277380.

The total contribution of these energies enabled the final estimation of the free binding energies in the complexes (Table 2). The lowest free binding energy contributing to more stability of the protein–ligand complex was observed by Marmin (−109.114 kJ/mol). Next amongst the four complexes was Pectachol (−63.487), followed by Karatavicinol (−57.644 kJ/mol), and Colladonin (−48.936 kJ/mol). The low binding energy of Marmin was much closer to that of CHEMBL1277380 (−111.732 kJ/mol) with that of Pectachol higher than ZINC8782981 (−54.399 kJ/mol). These energies address the potential of Marmin and Pectachol to bind most effectively at the active site of LmTR. LmTR–Marmin’s free binding energy correlated with the low binding energy (−9.3 kcal/mol) from docking. That of Pectachol showed a good free binding energy than that obtained from docking. This was better than that of Karatavicinol and Colladonin (Table 1).

3.15 Exploring possible implications and structure similarities of predicted leads

Karatavicinol and Marmin had lower binding energies of −9.4 and − 9.3 kcal/mol, respectively, as compared to Colladonin (−8.5 kcal/mol) and Pectachol (−8.5 kcal/mol). These binding energies are closer to that of FAD (−9.0 kcal/mol) for which can possibly compete in binding at the FAD domain. These compounds were concluded to have drug-likeness by satisfying Lipinski’s rule of 5. They also do not pass the blood–brain barrier which is good. Also, Marmin and Karatavicinol checked false for p-glycoprotein substrate. This gives the compounds an advantage to maintain their concentrations in cellular level to maximize efficacy. Pectachol and Colladonin however were implicated as P-gp substrates. These predicted preferable properties can favor their lead likeness and chances of going a long way in experimental studies. The four lead compounds were predicted as antileishmanial compounds. The four leads are confirmed not to be already existing antileishmanial drugs by structural similarity searches in but rather observed to be analogues of chrome 2-one. In regard to this, studies over the years have however shown some novel compounds such as 7-{[(2R*)-3,3-dimethyloxiran-2-yl]methoxy}-8-[(2R*,3R*)-3-isopropenyloxiran-2-yl]-2H-chromen-2-one and 7-methoxy-8-(4-methyl-3-furyl)-2H-chromen-2-one against Leishmania donovaniwith EC50 of 9.9 and 10.5 μg/mL, respectively [70]. These tested compounds with antileishmanial effect tend to be analogues of chromen-2-one. We emphasize that Karatavicinol is not a unique lead compound since it has already been experimented on other Leishmaniaspecies excluding L. major[57]. But the study identified it via these computational processes and therefore would report it as a potential compound against L. major. This augments the fact that the computational drug discovery pipeline has an optimistic potential of yielding good candidates for experimental work. Colladonin on the other hand is an enantiomer of Feselol for which Feselol is experimented as an antileishmanial agent [57]. Marmin also holds a very good potential of being an anti-ulcerative agent [71]. This favors it being a good compound for the treatment of cutaneous leishmaniasis. These compounds classified are coumarins and more other studies have reported good antileishmanial activities from this class of compounds [72, 73]. This work supports the fact that Karatavicinol, Marmin, Pectachol and Colladonin may possibly exhibit good antileishmanial activity if tested in vitro(Table A6).

Further in this study, the interaction of the active site residues with all four lead compounds showed hydrogen bonding with Val34, Thr51, Lys60, Thr160, Ala159, Arg287, Thr293, Asp327, Asn330, Thr335, and Gly376 (Table 1). Superimposition of the docked 2JK6 and co-crystallized revealed common residues such as Ser14, Gly15, Arg287, and Thr335 (Figure A1). These residues can be observed to be unique to the FAD domain of LmTR in anchoring the FAD molecule. Comparing these residues to the hydrogen bonding residues from the four leads shows that possible interruption of any of these residues can cause conformational changes which might not favor the selective binding of FAD at its domain. Baiocco et al. in 2009 identified Thr335 of trypanothione reductase at the FAD catalytic site of L. infantum[74]. They proposed that the FAD molecule binds tightly to the protein and orients itself towards the hydride transfer region of the active site by hydrogen bonding with specific residues Lys60, Thr335, and His461. Having observed this, interrupting one of these residues can potentially inhibit the reduction of T[S]2 by interfering with the hydride transfer. These compounds can potentially covey a competitive mode for binding to Thr355 which can affect the hydride transfer reaction in the active site preventing direct inactivation of trypanothione reductase. Other studies with quinone derivatives also have identified Thr355 and Ser14 as unique to the FAD domain of TR [75, 76].


4. Conclusion

Trypanothione reductase has been a well-investigated target essential for trypanosomatids. Its function in controlling oxidative stress in the parasite provided an opportunity to target the trypanothione biosynthesis pathway. A total of 11 hit compounds identified by pharmacophore modeling and virtual screening were filtered to four potential leads by considering their ADME with their molecular interactions in LmTR. MM-PBSA enabled the individual computation of active site residues that contributed significantly to binding. Efficient selective blockade of LmTR with these four coumarin compounds: Karatavicinol (7-[(2E,6E,10S)-10,11-dihydroxy-3,7,11-trimethyldodeca-2,6-dienoxy]chromen-2-one), Marmin (7-[(E,6R)-6,7-dihydroxy-3,7-dimethyloct-2-enoxy]chromen-2-one), Pectachol (7-[(6-hydroxy-5,5,8a-trimethyl-2-methylidene-3,4,4a,6,7,8-hexahydro-1H-naphthalen-1-yl)methoxy]-6,8-dimethoxychromen-2-one), and Colladonin (7-[[(4aS)-6-hydroxy-5,5,8a-trimethyl-2-methylidene-3,4,4a,6,7,8-hexahydro-1H-naphthalen-1-yl]methoxy]chromen-2-one) hold the potential to compromise the redox defenses of the parasites by inhibiting the FAD binding region and correspondingly increasing their sensitivity to redox-damage when carried out in in vitroand in vivostudies. Residues such as Asp35, Thr51, Lys61, Tyr198, and Asp327 are suspected to have critical role in the anchoring of FAD which contributes to the formation of reduced T[SH]2 in the reducing environment of amastigotes.



The authors are grateful to the West African Centre for Cell Biology of Infectious Pathogens (WACCBIP) at the University of Ghana for making Zuputo, a Dell EMC high-performance computing cluster, available for this study.


Conflict of interest

The authors declare no conflict of interest.


Leu10, Gly11, Gly13, Ser14, Gly15, Gly16, Val34, Asp35, Val36, Phe44, Ala46, Ala47, Gly50, Thr51, Cys52, Val55, Gly56, Cys57, Lys60, Lys61, Gly125, Phe126, Gly127, Ala128, Arg138, Ser140, Glu141, Pro143, Ala159, Thr160, Gly161, Ser162, Trp163, Pro164, Thr165, Thr177, Ser178, Asn179, Phe182, Tyr198, Ile199, Glu202, Phe203, Met282, Leu283, Ala284, Ile285, Gly286, Arg287, Arg290, Thr293, Leu294, Gln295, Ile325, Gly326, Asp327, Val332, Met333, Leu334, Thr335, Pro336, Val337, Ala338, Ile339, Asn340, Arg355, Thr357, Asp358, His359, Thr360, Lys361,Val362, Ala363, Cys364, Ala365, Phe367, Pro435, Glu436, Ile438, Gln439, Gly442, Ile443, Lys446

Table A1.

Predicted amino acid residues.

Figure A1.

(a) TheZ-score of the modeledLmTR (represented in dot) estimated to be −11.68 which is within the range of experimentally determined proteins by X ray method. (b) Verify 3D plot of the modeled protein structure,LmTR. This shows 91.65% of its amino acid residues with an average 3D-1D score greater than or equal to 0.2, which is a positive inference to the expected 80% of amino acids with 3D-1D score above or equal to 0.2.

Figure A2.

Ramachandran plot of the modeledLmTR protein structure. The percentage of residues in the most favored region (red) was 93.6% which is favorable for the protein’s stereochemistry. The percentage of residues in the allowed region (yellow) was 6%. Only 0.2% of protein residues (Phe45) showed probable stereochemical hindrance or collision.

Figure A3.

A 3D geometry of the generated pharmacophore. The nitrogen on the bicyclic ring of CNQB with the oxygen from the nitro group on its purine ring derivative contributed a hydrogen bond acceptor HBA (red sphere). The oxygen from the nitrogen dioxide group on the conjugated benzene in addition to the nitrogen on the five-member ring of PNTPC also contributed to HBA. Both had an aromatic ring (blue ring) as a common feature (blue ring in yellow 3D sphere) which contributed to hydrophobic interactions and the alkene feature shared amongst them generated the same hydrophobicity.

Figure A4.

(A) AUC score of 0.99 for the pharmacophore model. Determined at 1, 5, 10, and 100% of the selected database were the AUC and EF values as shown. The median is shown by dotted lines. If the curve is closer to the median it would suggest poor model. (B) AUC score of 0.702 generated for validating the docking system used. It verified the correlation between virtual screening performance and binding site descriptors of protein targets model (LmTR).

Figure A5.

(a) 2D schematic diagram of co-crystallized FAD (PDB ID: 2JK6) and FAD docked inLmTR superimposed together. Similar hydrogen bonding residues include Ser14, Gly15, Arg287, and Thr335. Similar hydrophobic residues in addition confirm the predicted active site. (b) Ligand alignment of co-crystalized FAD and FAD docked inLmTR.

NameP-fit scoreBinding energy (kcal/mol)Active/decoySource database
(S)-alpha-methoxy-alpha-trifluoromethyl-alpha-phenylacetate (MTPA)56.37−9.4ActiveNANPDB
4′-Methyl gossypetin56.17−8.2ActiveNANPDB
(+)-1,2-bis-(4-hydroxy-3-methoxyphenyl)-propane-1,3-diol [erythro form]55.9−7.4ActiveNANPDB
Drimartol A56.31−7.4ActiveNANPDB
3-(10-acetoxygeranyl)-4-acetoxy-p-coumaric acid56.14−7ActiveNANPDB

Table A2.

The 42 hits obtained from pharmacophore screening with their respective pharmacophore fit score, binding energies, and data sources.

Figure A6.

2D schematic diagram showing protein-ligands interaction of some leads at the active site ofLmTR. (A)LmTR–Karatavicinol interaction profile, (B)LmTR–Marmin interaction profile, (C)LmTR–Pectachol interaction profile, and (D)LmTR–Colladonin interaction profile.

Compound ZINC ID/nameNumber of Lipinski’s rules violatedMW (g/mol)No. HANo. HDxLogPWater solubility (mg/mL)Log SBio. Sc
ZINC954860810382.45524.52Moderately soluble−5.840.55
MTPA0470.52806.35Moderately soluble−6.000.55
Karatavicinol0400.51524.66Moderately soluble−4.850.55
Taccalin0418.4896−1.45Moderately soluble1.660.55
13-Hydroxyfeselol0400.51521.45Moderately soluble−5.930.55
Betaxanthin0370.4487−1.17Moderately soluble−2.110.55
Colladonin0384.51415.76Poorly soluble−6.500.55
Feselol0384.51415.76Poorly soluble−6.50.55
Pectachol0444.56615.70Poorly soluble−6.700.55

Table A3.

Physicochemical profiling of the 11 hit compounds.

All the hits showed good druglikeness. MW, molecular weight; No. HD, number of H-bond donors; Bio Sc, bioavailability score, No. HA, number of H-bond acceptors.

Compound ZINC IDGI absorptionBBB permeantP-gp substrateCYP1A2 inhibitorCYP2C19 inhibitorCYP2C9 inhibitorCYP2D6 inhibitorCYP3A4 inhibitor

Table A4.

Pharmacological profiling of top 11 compounds characterized by gastrointestinal (GI) absorption, blood brain barrier (BBB) permeant, p-gp substrates, and cytochrome P450 inhibitors.

Four compounds out of the 11 showed an appreciable pharmacological property. This included Karatavicinol, Marmin, Colladonin, and Pectachol.

Lead compoundsAntileishmanial predicted activity

Table A5.

This table shows the 10 top hit compounds and their predicted antileishmanial activity.

Karatavicinol, Taccalin, Marmin, 13-Hydroxyfeselol, Colladonin, Feselol, and Pectachol had a greater positive prediction above 0.5. If 0.5 < Pa < 0.7, the substance is likely to exhibit activity in experiment, but the probability of being a known pharmaceutical agent is less.

Figure A7.

Root mean square fluctuations of six complexes. The complexes are color coded in the graph. Karatavicinol and CHEMBL1277380 experienced the highest fluctuation at around residue number 80. Remaining complexes had similar patterns of fluctuations.

Figure A8.

The radius of gyration (Rg) plots of seven complexes within 100 ns simulation time. The complexes are represented in color code in the graph. Marmin showed the most preferentially well folded protein complex with Rg value of 2.33 nm.

Figure A9.

MM-PBSA plot of the binding free energy decomposition contribution per residue ofLmTR–Karatavicinol complex. Coded red lines represent surrounding active site amino acid residues.

Figure A10.

MM-PBSA plot of the binding free energy decomposition contribution per residue of LmTR–Pectachol complex. Coded red lines represent surrounding active site amino acid residues.

Figure A11.

MM-PBSA plot of the binding free energy decomposition contribution per residue ofLmTR–Colladonin complex. Coded red lines represent surrounding active site amino acid residues.

Figure A12.

Shows a graph of RMSD of the backbone of atoms in nm versus time in nanoseconds (ns). This graph is a representation of the average distance of the atoms of the residues at the backbone of the target protein. RMSD of 0.25 Å showed deviation from protein backbone.

Table A6.

Structure and IUPAC names of the three novel lead compounds.


A.1 List of abbreviations


absorption, distribution, metabolism, excretion and toxicity


area under curve


cytochromes P450


directory of useful (docking) decoys-enhanced


GROningen MAchine for Chemical Simulations


high performance computing



Log P

logarithm of the octan-1-ol/water partition coefficient


molecular dynamics


molecular mechanics Poisson Boltzmann surface area


molecular weight


permeability glycoprotein


prediction of activity spectra for substance


Protein Data Bank


radius of gyration


root mean square deviation


root mean square fluctuation


receiver operating characteristic


structure data file


simplified molecular input line entry system


universal force field


  1. 1. Murray HW, Berman JD, Davies CR, Saravia NG. Advances in leishmaniasis. Lancet. 2005;366(9496):1561-1577
  2. 2. Sundar S, Chakravarty J, Agarwal D, Rai M, Murray HW. Single-dose liposomal amphotericin B for visceral leishmaniasis in India. The New England Journal of Medicine. 2010;362(6):504-512
  3. 3. Ghisla S, Massey V. Mechanisms of flavoprotein-catalyzed reactions. European Journal of Biochemistry. 1989;181(1):1-17. DOI: 10.1111/j.1432-1033.1989.tb14688.x
  4. 4. Kinnings SL, Liu N, Tonge PJ, Jackson RM, Xie L, Bourne PE. A machine learning-based method to improve docking scoring functions and its application to drug repurposing. Journal of Chemical Information and Modeling. 2011;51(5):1195-1197. DOI: 10.1021/ci2001346
  5. 5. Shames SL, Fairlamb AH, Cerami A, Walsh CT. Purification and characterization of trypanothione reductase fromCrithidia fasciculata, a newly discovered member of the family of disulfide-containing flavoprotein reductases. Biochemistry. 1986;25(12):3519-3526
  6. 6. Krauth-Siegel RL, Enders B, Henderson GB, Fairlamb AH, Schirmer RH. Trypanothione reductase fromTrypanosoma cruzipurification and characterization of the crystalline enzyme. European Journal of Biochemistry. 1987;164(1):123-128
  7. 7. Khan MOF. Trypanothione reductase: A viable chemotherapeutic target for antitrypanosomal and antileishmanial drug design. Drug Target Insights. 2007;2:117739280700200
  8. 8. Beig M, Oellien F, Garoff L, Noack S, Krauth-Siegel RL, Selzer PM. Trypanothione reductase: A target protein for a combined in vitro and in silico screening approach. PLoS Neglected Tropical Diseases. 2015;9(6):e0003773. DOI: 10.1371/journal.pntd.0003773
  9. 9. Rodrigues RF, Castro-Pinto D, Echevarria A, Dos Reis CM, Del Cistia CN, Sant’Anna CMR, et al. Investigation of trypanothione reductase inhibitory activity by 1,3,4-thiadiazolium-2-aminide derivatives and molecular docking studies. Bioorganic and Medicinal Chemistry. 2012;20(5):1760-6
  10. 10. Veselovsky AV, Zharkova MS, Poroikov VV, Nicklaus MC. Computer-aided design and discovery of protein–protein interaction inhibitors as agents for anti-HIV therapy. SAR and QSAR in Environmental Research. 2014;25(6):457-471. DOI: 10.1080/1062936X.2014.898689
  11. 11. Taft CA, da Silva VB, da Silva CHT. Current topics in computer-aided drug design. Journal of Pharmaceutical Sciences. 2008;97(3):1089-1098
  12. 12. Ehrlich P. Über den jetzigen Stand der Chemotherapie. Berichte der Deutschen Chemischen Gesellschaft. 1909;42(1):17-47. DOI: 10.1002/cber.19090420105
  13. 13. Asimul Islam KK. Receptor chemoprint derived pharmacophore model for development of CAIX inhibitors. Journal of Carcinogenesis and Mutagenesis. 2014;s8(01):1-9
  14. 14. Chawla B, Madhubala R. Drug targets in leishmania. Journal of Parasitic Diseases. 2010;34(1):1-13. DOI: 10.1007/s12639-010-0006-3
  15. 15. Danesh-Bahreini MA, Shokri J, Samiei A, Kamali-Sarvestani E, Barzegar-Jalali M, Mohammadi-Samani S. Nanovaccine for leishmaniasis: Preparation of chitosan nanoparticles containing leishmania superoxide dismutase and evaluation of its immunogenicity in BALB/c mice. International Journal of Nanomedicine. 2011;6:835-842
  16. 16. Sundar S, Singh A. Chemotherapeutics of visceral leishmaniasis: Present and future developments. Parasitology. 2018;145(4):481-489
  17. 17. Croft SL, Seifert K, Yardley V. Current scenario of drug development for leishmaniasis. The Indian Journal of Medical Research. 2006;123(3):399-410
  18. 18. Sangshetti JN, Kalam Khan FA, Kulkarni AA, Arote R, Patil RH. Antileishmanial drug discovery: Comprehensive review of the last 10 years. RSC Advances. 2015;5(41):32376-32415
  19. 19. Agarwala R, Barrett T, Beck J, Benson DA, Bollin C, Bolton E, et al. Database resources of the National Center for Biotechnology Information. Nucleic Acids Research. 2018;46(D1):D8-D13
  20. 20. Pieper U, Webb BM, Dong GQ, Schneidman-Duhovny D, Fan H, Kim SJ, et al. ModBase, a database of annotated comparative protein structure models and associated resources. Nucleic Acids Research. 2014;42(D1):D336-D346. DOI: 10.1093/nar/gkt1144
  21. 21. Binkowski TA. CASTp: Computed atlas of surface topography of proteins. Nucleic Acids Research. 2003;31(13):3352-3355. DOI: 10.1093/nar/gkg512
  22. 22. Dundas J, Ouyang Z, Tseng J, Binkowski A, Turpaz Y, Liang J. CASTp: Computed atlas of surface topography of proteins with structural and topographical mapping of functionally annotated residues. Nucleic Acids Research. 2006;34:W116-W118
  23. 23. Dallakyan S, Olson AJ. Small-molecule library screening by docking with PyRx. Methods in Molecular Biology. 2015;1263:243-250. DOI: 10.1007/978-1-4939-2269-7_19
  24. 24. Padilha EC, Serafim RB, Sarmiento DYR, Santos CF, Santos CBR, Silva CHTP. New PPARα/γ/δ optimal activator rationally designed by computational methods. Journal of the Brazilian Chemical Society. 2016;27(9):1636-1647
  25. 25. DeLano WL. PyMOL: An open-source molecular graphics tool. CCP4 Newsletter on Protein Crystallography. 2002;40:82-92
  26. 26. Wiederstein M, Sippl MJ. ProSA-web: Interactive web service for the recognition of errors in three-dimensional structures of proteins. Nucleic Acids Research. 2007;35:W407-W410. DOI: 10.1093/nar/gkm290
  27. 27. Eisenberg D, Lüthy R, Bowie JU. [20] VERIFY3D: Assessment of protein models with three-dimensional profiles. Methods in Enzymology. 1997;277:396-404
  28. 28. Laskowski RA, MacArthur MW, Moss DS, Thornton JM. PROCHECK: A program to check the stereochemical quality of protein structures. Journal of Applied Crystallography. 1993;26(2):283-291
  29. 29. Patil R, Das S, Stanley A, Yadav L, Sudhakar A, Varma AK. Optimized hydrophobic interactions and hydrogen bonding at the target-ligand interface leads the pathways of drug-designing. PLoS One. 2010;5(8):e12029. DOI: 10.1371/journal.pone.0012029
  30. 30. Van Der Spoel D, Lindahl E, Hess B, Groenhof G, Mark AE, Berendsen HJC. GROMACS: Fast, flexible, and free. Journal of Computational Chemistry. 2005;26:1701-1718
  31. 31. Abraham MJ, Murtola T, Schulz R, Páll S, Smith JC, Hess B, et al. GROMACS: High performance molecular simulations through multi-level parallelism from laptops to supercomputers. SoftwareX. 2015;1–2:19-25
  32. 32. Wolber G, Langer T. LigandScout: 3-D pharmacophores derived from protein-bound ligands and their use as virtual screening filters. Journal of Chemical Information and Modeling. 2005;45(1):160-169. DOI: 10.1021/ci049885e
  33. 33. Bento AP, Gaulton A, Hersey A, Bellis LJ, Chambers J, Davies M, et al. The ChEMBL bioactivity database: An update. Nucleic Acids Research. 2014;42(D1):D1083-D1090. DOI: 10.1093/nar/gkt1031
  34. 34. Kim S, Thiessen PA, Bolton EE, Chen J, Fu G, Gindulyte A, et al. PubChem substance and compound databases. Nucleic Acids Research. 2016;44(D1):D1202-D1213. DOI: 10.1093/nar/gkv951
  35. 35. Irwin JJ, Shoichet BK. ZINC—A free database of commercially available compounds for virtual screening. Journal of Chemical Information and Modeling. 2005;45(1):177-182
  36. 36. Ntie-Kang F, Telukunta KK, Döring K, Simoben CV, Moumbock AFA, Malange YI, et al. NANPDB: A resource for natural products from Northern African sources. Journal of Natural Products. 2017;80(7):2067-2076. DOI: 10.1021/acs.jnatprod.7b00283
  37. 37. Mysinger MM, Carchia M, Irwin JJ, Shoichet BK. Directory of useful decoys, enhanced (DUD-E): Better ligands and decoys for better benchmarking. Journal of Medicinal Chemistry. 2012;55(14):6582-6594. DOI: 10.1021/jm300687e
  38. 38. Jiménez-Valverde A. Insights into the area under the receiver operating characteristic curve (AUC) as a discrimination measure in species distribution modelling. Global Ecology and Biogeography. 2012;21(4):498-507. DOI: 10.1111/j.1466-8238.2011.00683.x
  39. 39. O’Boyle NM, Banck M, James CA, Morley C, Vandermeersch T, Hutchison GR. Open babel: An open chemical toolbox. Journal of Cheminformatics. 2011;3(1):33. DOI: 10.1186/1758-2946-3-33
  40. 40. Goksuluk D, Korkmaz S, Zararsiz G, Karaagaoglu AE. easyROC: An interactive web-tool for ROC curve analysis using R language environment. The R Journal. 2016;8(2):213
  41. 41. Heifets A, Lilien RH. LigAlign: Flexible ligand-based active site alignment and analysis. Journal of Molecular Graphics & Modelling. 2010;29(1):93-101
  42. 42. Daina A, Michielin O, Zoete V. SwissADME: A free web tool to evaluate pharmacokinetics, drug-likeness and medicinal chemistry friendliness of small molecules. Scientific Reports. 2017;7(1):42717
  43. 43. Lipinski CA. Rule of five in 2015 and beyond: Target and ligand structural limitations, ligand chemistry structure and drug discovery project decisions. Advanced Drug Delivery Reviews. 2016;101:34-41
  44. 44. Lagunin A, Stepanchikova A, Filimonov D, Poroikov V. PASS: Prediction of activity spectra for biologically active substances. Bioinformatics. 2000;16(8):747-748. DOI: 10.1093/bioinformatics/16.8.747
  45. 45. Filimonov DA, Lagunin AA, Gloriozova TA, Rudik AV, Druzhilovskii DS, Pogodin PV, et al. Prediction of the biological activity spectra of organic compounds using the pass online web resource. Chemistry of Heterocyclic Compounds. 2014;50(3):444-457. DOI: 10.1007/s10593-014-1496-1
  46. 46. Kumari R, Kumar R, Lynn A. g_mmpbsa—A GROMACS tool for high-throughput MM-PBSA calculations. Journal of Chemical Information and Modeling. 2014;54(7):1951-1962. DOI: 10.1021/ci500020m
  47. 47. Uziela K, Shu N, Wallner B, Elofsson A. ProQ3: Improved model quality assessments using Rosetta energy terms. Scientific Reports. 2016;6(1):33509
  48. 48. Shamsara J. Correlation between virtual screening performance and binding site descriptors of protein targets. International Journal of Medicinal Chemistry. 2018;2018:1-10
  49. 49. Fu Y, Zhao J, Chen Z. Insights into the molecular mechanisms of protein-ligand interactions by molecular docking and molecular dynamics simulation: A case of oligopeptide binding protein. Computational and Mathematical Methods in Medicine. 2018;2018:1-12
  50. 50. Bhakat S. Effect of T68A/N126Y mutations on the conformational and ligand binding landscape of Coxsackievirus B3 3C protease. Molecular BioSystems. 2015;11(8):2303-2311
  51. 51. Barthe L, Woodley J, Houin G. Gastrointestinal absorption of drugs: Methods and studies. Fundamental & Clinical Pharmacology. 1999;13(2):154-168. DOI: 10.1111/j.1472-8206.1999.tb00334.x
  52. 52. Smith HS. Opioid metabolism. Mayo Clinic Proceedings. 2009;84(7):613-624
  53. 53. Lespine A, Ménez C, Bourguinat C, Prichard RK. P-glycoproteins and other multidrug resistance transporters in the pharmacology of anthelmintics: Prospects for reversing transport-dependent anthelmintic resistance. International Journal for Parasitology: Drugs and Drug Resistance. 2012;2:58-75
  54. 54. Pardridge WM. Drug transport across the blood–brain barrier. Journal of Cerebral Blood Flow and Metabolism. 2012;32(11):1959-1972. DOI: 10.1038/jcbfm.2012.126
  55. 55. Lagunin A, Filimonov D, Poroikov V. Multi-targeted natural products evaluation based on biological activity prediction with PASS. Current Pharmaceutical Design. 2010;16(15):1703-1717
  56. 56. Filimonov DA, Druzhilovskiy DS, Lagunin AA, Gloriozova TA, Rudik AV, Dmitriev AV, et al. Computer-aided prediction of biological activity spectra for chemical compounds: Opportunities and limitation. Biomedical Chemistry: Research and Methods. 2018;1(1):e00004
  57. 57. Amin A, Tuenter E, Cos P, Maes L, Exarchou V, Apers S, et al. Antiprotozoal and antiglycation activities of sesquiterpene coumarins from ferula narthex exudate. Molecules. 2016;21(10):1287
  58. 58. Gliszczyńska A, Brodelius PE. Sesquiterpene coumarins. Phytochemistry Reviews. 2012;11(1):77-96. DOI: 10.1007/s11101-011-9220-6
  59. 59. Karplus M, McCammon JA. Molecular dynamics simulations of biomolecules. Nature Structural Biology. 2002;9(9):646-652. DOI: 10.1038/nsb0902-646
  60. 60. Hollingsworth SA, Dror RO. Molecular dynamics simulation for all. Neuron. 2018;99(6):1129-1143
  61. 61. Cheng X, Ivanov I. Molecular dynamics. Methods in Molecular Biology. 2012;929:243-285
  62. 62. Sinha S, Wang SM. Classification of VUS and unclassified variants in BRCA1 BRCT repeats by molecular dynamics simulation. Computational and Structural Biotechnology Journal. 2020;18:723-736
  63. 63. Perez A, Morrone JA, Simmerling C, Dill KA. Advances in free-energy-based simulations of protein folding and ligand binding. Current Opinion in Structural Biology. 2016;36:25-31
  64. 64. Kwofie SK, Dankwa B, Enninful KS, Adobor C, Broni E, Ntiamoah A, et al. Molecular docking and dynamics simulation studies predict Munc18b as a target of mycolactone: A plausible mechanism for granule exocytosis impairment in buruli ulcer pathogenesis. Toxins. 2019;11(3)
  65. 65. Sartori GR, Nascimento AS. Comparative analysis of electrostatic models for ligand docking. Frontiers in Molecular Biosciences. 2019;6(52). DOI: 10.3389/fmolb.2019.00052/full
  66. 66. Deng N, Zhang P, Cieplak P, Lai L. Elucidating the energetics of entropically driven protein–ligand association: Calculations of absolute binding free energy and entropy. The Journal of Physical Chemistry. B. 2011;115(41):11902-11910. DOI: 10.1021/jp204047b
  67. 67. Fiorentini R, Kremer K, Potestio R. Ligand-protein interactions in lysozyme investigated through a dual-resolution model. Proteins: Structure, Function, and Bioinformatics. 2020;88(10):1351-1360. DOI: 10.1002/prot.25954
  68. 68. Campanera JM, Pouplana R. MMPBSA decomposition of the binding energy throughout a molecular dynamics simulation of amyloid-beta (Aß10−35) aggregation. Molecules. 2010;15(4):2730-2748
  69. 69. Genheden S, Ryde U. The MM/PBSA and MM/GBSA methods to estimate ligand-binding affinities. Expert Opinion on Drug Discovery. 2015;10(5):449-461. DOI: 10.1517/17460441.2015.1032936
  70. 70. Arango V, Robledo S, Séon-Méniel B, Figadère B, Cardona W, Sáez J, et al. Coumarins fromGalipea panamensisand their activity againstLeishmania panamensis. Journal of Natural Products. 2010;73(5):1012-1014. DOI: 10.1021/np100146y
  71. 71. Alam F. Anti-ulcer plants from North-East India—A review. Der Pharmacia Lettre. 2019;11(6):73-96
  72. 72. Sajjadi S, Eskandarian A-A, Shokoohinia Y, Yousefi H-A, Mansourian M, Asgarian-Nasab H, et al. Antileishmanial activity of prenylated coumarins isolated fromFerulago angulataandPrangos asperula. Research in Pharmaceutical Sciences. 2016;11(4):324
  73. 73. Ferreira ME, Rojas de Arias A, Yaluff G, de Bilbao NV, Nakayama H, Torres S, et al. Antileishmanial activity of furoquinolines and coumarins fromHelietta apiculata. Phytomedicine. 2010;17(5):375-378
  74. 74. Baiocco P, Colotti G, Franceschini S, Ilari A. Molecular basis of antimony treatment in leishmaniasis. Journal of Medicinal Chemistry. 2009;52(8):2603-2612. DOI: 10.1021/jm900185q
  75. 75. Ravi Kumar G, Jagannadham Medicherla V. Molecular docking based inhibition of trypanothione reductase activity by taxifolin novel target for antileishmanial activity. Journal of Applied Pharmaceutical Science. 2012;2(10):133-136
  76. 76. Venkatesan SK, Shukla AK, Dubey VK. Molecular docking studies of selected tricyclic and quinone derivatives on trypanothione reductase ofLeishmania infantum. Journal of Computational Chemistry. 2010;31(13):2463-2475. DOI: 10.1002/jcc.21538

Written By

Samuel K. Kwofie, Gabriel B. Kwarko, Emmanuel Broni, Michael B. Adinortey and Michael D. Wilson

Submitted: August 10th, 2021 Reviewed: September 23rd, 2021 Published: December 18th, 2021