|Home | About | Journals | Submit | Contact Us | Français|
Plasmodium falciparum dihydroorotate dehydrogenase (PfDHODH), a key enzyme in the de novo pyrimidine biosynthesis pathway, which the Plasmodium falciparum relies on exclusively for survival, has emerged as a promising target for antimalarial drugs. In an effort to discover new and potent PfDHODH inhibitors, 3D-QSAR pharmacophore models were developed based on the structures of known PfDHODH inhibitors and the validated Hypo1 model was used as a 3D search query for virtual screening of the National Cancer Institute database. The virtual hit compounds were further filtered based on molecular docking and Molecular Mechanics/Generalized Born Surface Area binding energy calculations. The combination of the pharmacophore and structure-based virtual screening resulted in the identification of nine new compounds that showed >25% inhibition of PfDHODH at a concentration of 10 μM, three of which exhibited IC50 values in the range of 0.38–20 μM. The most active compound, NSC336047, displayed species-selectivity for PfDHODH over human DHODH and inhibited parasite growth with an IC50 of 26 μM. In addition to this, thirteen compounds inhibited parasite growth with IC50 values of ≤ 50 μM, four of which showed IC50 values in the range of 5–12 μM. These compounds could be further explored in the identification and development of more potent PfDHODH and parasite growth inhibitors.
Malaria is a life-threatening parasitic disease, which continues to cause a severe global health problem. The World Health Organization estimated that in 2013 alone, malaria infected 198 million people and resulted in over half a million deaths.1 Malaria mostly affects children living in Africa, where it is believed that a child dies from malaria every minute.1 The disease is caused by Plasmodium parasites, which are spread to people through the bites of infected Anopheles mosquitoes. Of the five Plasmodium parasite species that cause malaria in humans, Plasmodium falciparum (P. falciparum) is responsible for most of the severe clinical malaria cases. No effective vaccines are available to date and thus antimalarial control programs and chemotherapy are used for prevention and treatment of this disease, respectively.2 Unfortunately, the chemotherapy-based treatments using traditional antimalarial drugs face widespread drug-resistance, compromising their efficacy. For example, some commonly used antimalarial drugs such as chloroquine 3, atovaquone 4, pyrimethamine 5 and sulfadoxine 6 have been withdrawn in many areas due to parasite resistance. Moreover, resistance has been reported to artemisinin-based combination therapies (ACTs), a new treatment option to combat drug-resistance, in the Thai-Cambodian border. 7 This highlights the pressing need for the development of novel, non-cross-resistant and effective antimalarial drugs. The identification of unique biochemical processes that are critical for parasite survival is an important step in the identification of novel drug targets for new antimalarial drugs.8
Pyrimidines are required for many biochemical processes including DNA and RNA synthesis, protein glycosylation and membrane lipid synthesis. 9 These essential precursor molecules are synthesized in many organisms, including humans, by both de novo pyrimidine biosynthetic pathways as well as salvage pathways that recover purine and pyrimidine bases formed during nucleic acid degradation. In contrast, the P. falciparum parasite genome lacks the required components for the pyrimidine salvage pathway 10 and thus the parasite relies exclusively on the de novo pyrimidine biosynthetic pathway.11 For this reason, the de novo pyrimidine biosynthetic pathway of P. falciparum has become an attractive target for the development of novel therapeutics for malaria. DiHydroOrotate DeHydrogenase (DHODH), the fourth key enzyme in P. falciparum de novo pyrimidine biosynthesis, catalyzes the oxidation of dihydroorotate to produce orotate in the presence of the co-factors flavin mononucleotide (FMN) and ubiquinone (CoQ). 12-14 P. falciparum DHODH (PfDHODH) and human DHODH (hDHODH) are members of the type II family of DHODH and are localized in the outer side of the inner mitochondrial membrane. They use respiratory quinones exclusively as terminal electron acceptors. 13, 15, 16 Recent studies showed that the physiological oxidant ubiquinone, provided by the mitochondrial electron transport chain, is needed by PfDHODH for pyrimidine biosynthesis. Studies have also shown that PfDHODH inhibition leads to parasite death in both cell culture and animal models. 17, 18 Furthermore, inhibitors of hDHODH are presently in use for the treatment of rheumatoid arthritis 19 and have been evaluated as anti-tumor 20 and immunosuppressive agents 21, demonstrating that DHODH is a druggable target.
Crystal structures of PfDHODH in complex with a number of diverse inhibitors have been determined 17, 22-27, providing insights into the structural basis for inhibition. PfDHODH has a large flexible active site that can accommodate structurally diverse inhibitors 26. Key hydrogen-bond interactions between the inhibitor and His185 and Arg265 residues are important for inhibitor binding and stacking interactions between Phe227 and Phe88 residues leads to favorable binding of inhibitors with large aromatic rings. Deng et al 25, 26 elucidated the structural basis for the species selective binding of the triazolopyrimidines to PfDHODH, identifying key amino acid residues that confer selectivity. Analysis of hDHODH and PfDHODH crystal structures revealed that replacement of Ala59 and Pro364 residues in hDHODH with the more bulky Phe188 and Met536 residues in PfDHODH hinders the hydrophobic pocket that binds the biphenyl moiety of brequinar in hDHODH. In addition, the substitution of Thr63 and Met111 residues in hDHODH with Gly192 and Leu240 residues in PfDHODH leads to the formation of the pocket that binds the naphthyl group of triazolopyrimidines in PfDHODH.28
In the past decade, several classes of PfDHODH inhibitors with diverse chemical structures have been developed (Figure 1). Interestingly, modification of known hDHODH inhibitors using de novo design led to the identification of a new PfDHODH inhibitor (compound 1, Figure 1) 29. Phillips and Rathod et al. identified numerous classes of potent and selective PfDHODH inhibitors through PfDHODH enzyme-based screening 18, 30 and a lead optimization program on the selected triazolopyrimidine-based compound 2 led to the discovery of DSM265 (compound 3) which is currently in clinical development. 23, 31-33 In addition, Clardy and Wirth et al. identified some selective PfDHODH inhibitors using target-based high throughput screening with the Genzyme library; compound 4 was the most potent among them. 34 Xu et al. discovered the hit compound 5 by structure-based screening with the SPECS database and subsequent lead optimization led to the potent and selective compound 6. 35 These studies demonstrate that the design of potent species-selective PfDHODH inhibitors is feasible using a combination of biochemical, structural and computational approaches.
Diverse computational tools are systematically integrated to speed up and facilitate hit identification, hit-to-lead selection, and lead optimization in drug discovery process.36-38 In such an integration of computational methods to identify new hits, pharmacophore modelling is initially employed as a powerful search tool to retrieve active compounds from databases. This search retrieves compounds that are similar to the pharmacophore model or entirely diverse scaffolds as the features of a pharmacophore model can map to multiple structural elements of the functional groups of compounds. The high scored compounds are collected as hits as they better map to the pharmacophore model.39, 40 This step significantly decreases the number of hits. As the second step, the best scored compounds retrieved by pharmacophore-based screening are subsequently docked into the binding pocket of a target to narrow the number of active hits.41, 42 The combination of pharmacophore- and docking-based virtual screening has been shown to be beneficial as it includes all possible information.43-45
In this study, 3D quantitative structure–activity relationship (3D-QSAR) pharmacophore models were built based on a series of dihydrothiophenone derivatives that show PfDHODH inhibition. The validated model was used as a 3D search query for screening the National Cancer Institute (NCI) compound database to identify new inhibitors. Subsequently, the hits with good pharmacophoric fit values were docked into the crystal structure of PfDHODH using Glide standard precision (SP) and extra precision (XP) docking algorithms and filtered accordingly. The hit compounds were further evaluated based on their predicted binding affinity for PfDHODH calculated using Prime/MM-GBSA. Following this, sixty two hits were finally selected for experimental validation. Nine hit compounds exhibited >25 % inhibition against PfDHODH, with three of these showing IC50 values between 0.38 and 20 μM. A hDHODH inhibition experiment confirmed the species-selectivity of the most active compound. In addition to this, thirteen compounds inhibited the growth of P. falciparum NF54 strain with IC50 values of ≤ 50 μM, with four of them showing IC50 values between 5–12 μM.
A set of 38 PfDHODH inhibitors were collected from recently published literature reported by Xu et al. 35 The inhibitors were derived from the dihydrothiophenone scaffold and the in vitro bioactivities of the collected inhibitors were expressed as the concentration of the test compounds that inhibited the activity of PfDHODH by 50% (IC50). The two-dimensional (2D) chemical structures of the inhibitors were sketched using ChemDraw Ultra and saved in MDL (Molecular Design Limited) mol file format. Subsequently, they were imported into Discovery Studio 4.0 (DS4) (Dassault Systèmes BIOVIA, Discovery Studio Modeling Environment, Release 4.0, San Diego: Dassault Systèmes, 2015) and converted into corresponding standard 3D-structures, which were utilized as starting conformations for conformation generations.
According to the Catalyst program guidelines, the initial group of 38 inhibitors was then divided into a training set (Figure S1) and a test set (Figure S2). The training set consisted of 19 structurally diverse compounds with bioactivities ranging from 6 to 39,450 nM, spanning over four orders of magnitude, including several most active, moderately active, and least active compounds. The activity values were classified as follows: IC50 ≤ 1000 nM means the compounds are the most active (represented as +++); 1000 nM < IC50 ≤ 10,000 nM means the compounds are the moderately active (represented as ++); and IC50 > 10,000 nM means the compounds are the least active (represented as +). The training set was then used to generate 3DQSAR pharmacophore models. The remaining 19 compounds were used as a test set to evaluate the predictive capability of the pharmacophore models.
The molecular flexibility of each compound in the training and test sets was modelled by generating multiple conformers within a specific energy range. For this, 250 conformational models for each compound were generated using the Poling algorithm's “best quality” conformational search option 46 within the ConFirm module of DS4, which uses generalized CHARMm force fields. 47 An energy threshold of 15 kcal/mol was used to ensure maximum coverage of the conformational space; default settings were used for all other parameters.
3D-QSAR pharmacophore hypotheses were generated using the HypoGen algorithm of DS4. In the HypoGen pharmacophore model generation, a minimum of 0 to a maximum of 5 features were selected and used to build a series of hypotheses using an uncertainty value of 2. HypoGen generates quantitative pharmacophore models based on the most active compounds from the training set and examines the training set compounds to correlate the structure-activity relations. In addition, the activity of each compound in the training set is estimated using the regression parameter, which is calculated by plotting the geometric fit value versus the negative logarithm of activity for each compound. The top 10 hypotheses with significant statistical parameters were selected on the basis of a high correlation coefficient (r), low total cost, and low root mean square (RMS).
Each pharmacophore model's (Hypo1-Hypo10) capability of predicting the activity of external compounds was determined using the test set compounds. The test set compounds were prepared in the same way as the training set compounds and the performance of the models was examined by regressing them against the test set compounds.
The statistical significance of a pharmacophore model can be described in terms of the fixed cost, null cost and total cost. 48 The fixed cost is the simplest model that fits all data perfectly and the null cost is the highest cost of a pharmacophore with no features and estimates activity to be the average of the activity data of the training set compounds. The total cost of any pharmacophore hypothesis should be close to the fixed cost to provide valuable models. In addition, the pharmacophore model is considered significant when the difference between the null and fixed cost value is large. A cost difference value should be 40–60 bits for a pharmacophore hypothesis, indicating that the model has 75–90% probability of correlating the data. The configuration cost and error cost, determine the quality of any pharmacophore hypotheses with predictive values. The configuration cost represents the complexity of the pharmacophore hypothesis space and should have a value of less than 17. The error cost dependents on the root-mean-square deviations (RMSDs) between the estimated and the experimental activities of the training set molecules. The RMSD represents the quality of the correlation between the experimental and the predicted activity data. The best pharmacophore model should have the highest cost difference, the lowest RMSD, and the best correlation coefficient.
The cross-validation was performed using the Fischer's randomization test 49 which uses the CatScramble program to randomize the experimental activities of the training set. The confidence level was set to 95% and 19 different random spreadsheets were generated and used to construct hypotheses using exactly the same features and parameters used in generating the original pharmacophore hypotheses. This was done to ensure that the cost values of the randomized data set were not similar or better than the cost values of the original hypotheses.
The validated pharmacophore model, Hypo1, containing four pharmacophore features, was used as a 3D-search query to retrieve lead-like compounds from the NCI database containing 265,242 compounds (released in May, 2012, http://cactus.nci.nih.gov/download/nci/). The database was first prepared and duplicates were removed using the Prepare Ligands protocol of DS4, which returned 246,477 compounds. The molecular flexibility of each compound in the database was then modelled by generating multiple conformers with “best quality” conformational search option. The best/flexible search option of Catalyst program was applied to pharmacophore-based virtual screening to retrieve compounds from the database. A molecule was only retrieved as a hit if it fitted to all the features of a pharmacophore model. The hit compounds were ranked according to the fit value and the compounds with good fit values were docked into the crystal structure of PfDHODH.
The crystal structure of PfDHODH in complex with inhibitor DSM1 (PDB: 3I65) was obtained from the protein data bank (PDB) and prepared for structure-based virtual screening 26. All HETATM residues and crystal water molecules were deleted except a crucial water molecule (W15), which mediates a hydrogen bond interaction between DSM1 and the binding site. The protein was prepared using Protein Preparation Wizard in Schrödinger 9.7 (Schrödinger, Inc., New York, NY, USA). This wizard was used to correct bond orders, add hydrogen atoms, create zero-order bonds to metals, optimize the orientations of added hydrogen for optimal hydrogen bond formation, and finally to minimize heavy atoms to a RMSD threshold of 0.3 Å using OPLS_2005 (optimized potential for liquid simulations_2005) force fields.
The Glide docking program in Schrodinger 9.7 (Schrödinger, Inc., New York, NY, USA) was used for the docking experiments. The docking method described below was validated by re-docking the DMS1 structure into the PfDHODH crystal structure and calculating the RMSD between the top docked pose and the bound DMS1 conformation in the crystal structure. A docking grid box on the centroid of the DSM1 in the crystal structure was generated using Glide-Receptor Grid Generation with default parameters for van der Waals (VDW) radius scaling. The hit molecules were optimized using the LigPrep module of the Schrödinger suite with OPLS_2005 force fields. This module generates possible 3D conformations for each ligand with various ionization states at pH 7.0 ± 2.0, tautomers, stereochemistries and ring conformations. All the generated conformations were saved as maestro files and used for docking. Standard precision (SP) Glide docking (default parameters) was used to dock hit molecules into the DSM1 binding site and the best pose for each hit was chosen based on the Glide docking score (G-Score). The hits with the best poses were then subjected to Extra precision (XP) Glide docking, a more precise scoring function, and finally, the top ranked poses were retained for Molecular Mechanics/Generalized Born Surface Area (MM-GBSA) calculation.
The binding free energies of docking poses, obtained from XP Glide docking, were calculated using the MM-GBSA method (Prime version 2.1, Schrödinger, LLC, New York, NY, 2009) with default setting. This method calculates energies using the OPLS-AA force fields for molecular mechanics energy (EMM), the surface-generalized born model for polar solvation energy (GSGB), and a nonpolar solvation term (GSA). The binding free energy (ΔGbind) was calculated as follows:
The scored poses by Prime MM-GBSA were visually inspected for protein-inhibitor interactions similar to those seen in inhibitor bound PfDHODH crystal structures.
The compounds that showed favorable interactions with PfDHODH were compared with eight known PfDHODH inhibitors for structural similarity (Figure S3) using SciTegic extended-connectivity fingerprints (ECFP_4) in DS4. The ECFP_4 is an atom type-based method, which computes the fingerprint features for each atom that is within a diameter of four bonds to the neighbor atoms. Tanimoto (Tc) similarity distance was used as a metric to filter out compounds that were structurally similar to the known PfDHODH inhibitors. On the basis of the metric, compounds were selected and clustered using Cluster Ligands protocol in DS4. After examining the favorable interaction and structural diversity, 39 compounds were finally obtained from the NCI and tested in vitro on PfDHODH and PfNF54.
The PfDHODH-NSC336047 structure was superimposed with the docked hDHODH-NSC336047, PfDHODH-DSM1 (PDB: 3I65)26, PfDHODH-A771726 (PDB 1TV5)50, hDHODH-Brequinar (PDB 1D3G)51 structures by aligning the backbone atoms of full-length structures in PyMOL.52 Structures were displayed using the graphics program PyMOL.
The structural similarity of identified hit compounds was analyzed using SciFinder (SciFinder; American Chemical Society, 2015) and compared to known PfDHODH inhibitors found in ChEMBL53 and 100 known active PfDHODH inhibitors (Table S2). Tanimoto similarity of each hit compound was calculated against ChEMBL and known active PfDHODH inhibitors using the “Find similar molecules by fingerprints” protocol in DS4, along with use of the ECFP_4 fingerprint. This protocol allowed us to calculate the minimum, maximum, and average similarities as measures of proximity to known PfDHODH inhibitors.
P. falciparum and human DHODH were expressed as recombinant proteins in BL21 (DE3) phage-resistant cells and purified as previously described.18, 30 Steady-state kinetic analysis for IC50 determination was performed using the colorimetric assay that monitors the reduction of 2, 6-dichloroindophenol (DCIP) at 600 nm (ε = 18.8 mM) for measuring DHODH inhibition. The assay was carried out using a buffer containing 100 mM HEPES, pH 8.0, 150 mM NaCl, 10 % Glycerol, 0.1 % Triton X-100, 20 μM CoQD (coenzyme QD), 200 μM L-dihydroorotate and 120 μM DCIP. The enzymatic reaction was initiated by the addition of 5 nM of enzyme. Data were fitted to the log [I] vs response (three parameters) equation (Y=Bottom + (Top-Bottom)/(1+10^((X-LogIC50))) or to the standard IC50 equation (Y=1/(1+X/(IC50))) for compounds with IC50 >10 μM in Graph Pad Prism.
The in vitro antiplasmodial activity screening was performed against chloroquine sensitive (CQS) P. falciparum NF54 strain at the University of Cape Town (Division of Clinical Pharmacology) for selected thirty nine virtual hit compounds. The screening method was based on the parasite lactate dehydrogenase assay using chloroquine and artesunate as reference standards according to previously described methods54, 55 (Supporting information). An additional twenty three compounds were tested at the Swiss Tropical and Public Health (Swiss TPH) Institute. The screening method was based on the [3H] hypoxanthine incorporation assay using chloroquine and artesunate as the standard reference drugs according to previously described procedures 56, 57 (Supporting information).
38 PfDHODH inhibitors were divided into a training set of compounds (Figure S1) and a test set of compounds (Figure S2) for 3D-QSAR pharmacophore model generation and validation, respectively. The training set compounds play a key role in determining the quality of the pharmacophore models generated; while the test set compounds served to evaluate the predictive capability of the resultant pharmacophore models. Both the training set and the test set were made up of molecules with a large range of activities ensuring that critical information on the pharmacophoric requirements for PfDHODH inhibition could be incorporated into the model and evaluated.
The training set, consisting of 19 structurally diverse dihydrothiophenone inhibitors with inhibitory activities spanning over four orders of magnitude (IC50 6–39,450 nM; Figure S1), was used to generate pharmacophore models for PfDHODH inhibitors. The cost values, correlation coefficients (R), root-mean square (RMS) values and the pharmacophore features of the top 10 ranked hypotheses (Hypo1–Hypo10) are shown in Table 1.
The cost analysis statistical calculation was used to determine the quality of the pharmacophore hypotheses. The null cost, fixed cost, and configuration costs for the top 10 scored hypotheses were 198.74, 68.27 and 11.99, respectively. The magnitude of the cost difference between the null cost (198.74) and the total cost (87.76) for the top ranked hypothesis, Hypo1, was 110.98; a cost difference greater than 60 indicates that the model has good predictive power and should represent 90% of the true correlation. In addition, the total cost (87.76) of the Hypo1 is close to the fixed cost (68.27), confirming that the generated Hypo1 was not obtained by chance. The configuration cost value must be below 17 for the pharmacophore model to be considered as a good hypothesis. In our case a configuration cost value of 11.991 was obtained for all hypotheses (Table 1), indicating that the generated pharmacophore hypotheses are reasonable.
All 10 of the hypotheses were evaluated using both the training and the test set compounds. The Hypo1 model gave the lowest RMS (1.385) and the best correlation between the predicted activities and the experimental activities for the test set compounds. The evaluation of the Hypo1 model is described in detail below. Hypo1 consists of four features: two hydrogen-bond acceptors (A), one hydrophobic-aliphatic (H) feature and one ring-aromatic (R) feature. All 10 hypotheses had the same pharmacophore features as Hypo1. However, their spatial locations and statistical parameters differed (Table 1 and Figure S4). The 3D-space and distance constraints of the pharmacophore features of Hypo1 are shown in Figure 2.
The activities predicted by the Hypo1 pharmacophore model for the 19 training set compounds, together with their experimental activities, are listed in Table 2. These compounds were classified into three groups (activity scales) based on their experimental activities (IC50): most active (IC50 < 1000 nM, +++); moderately active (1000 nM ≤ IC50 < 10000 nM, ++); least active (IC50 ≥ 10000 nM, +). The activities predicted by Hypo1 for all the compounds fell within the correct activity scales. The ratio between experimental and predicted activities (error) for most of the compounds was less than 5 (Table 2) and the correlation between the experimental and predicted activities was 0.935 (Figure S5), signifying good consistency between predicted and experimental IC50 values.
To further evaluate the robustness of the Hypo1 model, the training set compounds were mapped onto the Hypo1 model using the best fit option; the pharmacophoric fit values for each compound are reported in Table 2. The most active compound (compound 7, IC50=6 nM) mapped well to all the pharmacophore features of the model [Figure 3(A)] with a good pharmacophoric fit value of 9.63 and a predicted IC50 of 20.47 nM. On the contrary, while the least active compound (compound 25, IC50=39,450 nM) mapped onto the ring-aromatic and hydrophobic-aliphatic features, it only mapped onto one of the two hydrogen-bond acceptor features [Figure 3(B)] and thus exhibited a lower pharmacophoric fit value of 6.61. The predicted IC50 value for this compound was over 21,685.3 nM, 3 orders of magnitude higher than the IC50 features predicted for the most active compound (Table 2), suggesting that a second hydrogen-bond acceptor is important for PfDHODH binding. However, it cannot be ruled out that the presence of an intramolecular hydrogen-bond in 25, as suggested by Xu et al 35, may control its conformation to miss one of the two hydrogen-bond acceptor features.
The Hypo1 model's capability to predict the activity of inhibitors that are not part of the training set was evaluated using the 19 PfDHODH inhibitors in the test set (Figure S2). The IC50 values predicted by the Hypo1 model for the test set compounds together with their experimental IC50 values and the error values for the predictions are shown in Table S1. Of the 19 test set compounds, 17 had error values below 10; the other two compounds (compounds 31 and 32) had error values of 15.90 and 11.14, respectively. There was a good correlation (R = 0.933) between the experimental activities and the predicted activities for the test set (Figure S5). These results demonstrate that the Hypo1 model is capable of predicting activities of PfDHODH inhibitors not included in the training set.
The Fischer's randomization test was used to evaluate the statistical significance of the Hypo1 model. This cross-validation method produced 19 random pharmacophore spreadsheets, at 95% confidence level, using the same features and parameters as used in the generation of the original 10 pharmacophore hypotheses. The cost value of Hypo1 is lower than the 19 randomly generated hypotheses (Figure 4), indicating that the Hypo1 model is superior and was not generated by chance. Furthermore, the result of Fischer's test confirms that there is a 95% chance that the Hypo1 represents a true correlation for the training set compounds, providing further evidence that the pharmacophore model, Hypo1, is a reliable model with statistical significance.
While this manuscript was in preparation, Aher and Roy reported a pharmacophore model generated using a similar series of dihydrothiophenone PfDHODH inhibitors. 58 Their model differed from the Hypo1 model presented here both in terms of features as well as feature constraints, which may be due to the different training set compounds used. 59 In contrast to their model, our model has an additional hydrophobic-aliphatic feature that maps to the ethyl group of many of the active compounds as shown in Figure 3. In addition, contrary to their validation, we have adopted the standard regression method on the test set, which showed good correlation (Figure S5). Furthermore, we have experimentally demonstrated in the following section that our model is capable of picking the active molecules from the NCI database, suggesting that the model is reasonable and more reliable.
Virtual screening, a complementary approach to the experimental high-throughput screening, has been demonstrated to increase the success rate in the hit identification stage of the drug discovery process.60 Many studies have suggested that combining pharmacophore-based screening with structure-based methods could enhance the enrichment of virtual screening. In this study, a multistep virtual screening protocol, combining pharmacophore-based screening with structure-based methods, was adopted to identify new and potent PfDHODH inhibitors in the NCI database. A schematic representation of the screening cascade is shown in Figure 5.
The NCI open library has previously been used to find new anti-cancer and anti-viral agents using both in vitro and in vivo screening methods. 61 In this study, the NCI database, consisting of 246,477 unique compounds, was screened for PfDHODH inhibitors using the pharmacophore model Hypo 1 as a 3D-search query. The database was spiked with 100 known active PfDHODH inhibitors (≤ 1μM) (Table S2), which were not used to generate the Hypo1 model, and a multi-conformer database was generated. The model retrieved 50,211 hit compounds from the 246,577 compounds, of which 95 were known inhibitors. The enrichment factor (EF = [95 × (246577+100)] ÷ (50211 × 100) = 4.66) was calculated to be 5, indicating that the Hypo 1 model was 5 times more likely to pick active molecules from the database than would be expected by chance. 62 Table 3 lists the statistical parameters from the pharmacophore-based virtual screening. To increase the stringency of the search, a pharmacophore fit value cut-off of ≥ 9, based on the pharmacophore fit of the most active compound in the training set, was applied. Therefore, the compounds with a fit value of <9 were considered as inactive and were omitted. This criterion reduced the number of active hits to 3529. To reduce false positives and to obtain lead-like compounds, these hits were subjected to docking-based virtual screening.
Initial docking-based virtual screening was performed using the Glide SP method, which enabled the rapid docking of the 3529 hit compounds, obtained from the pharmacophore-based screening, into the PfDHODH crystal structure (PDB ID: 3I65). The capability of the docking protocol was validated prior to screening the hit molecules by re-docking the DSM1 into the PfDHODH crystal structure to see whether or not the crystal bound conformation of the ligand could be reproduced. The heavy atom RMSD between the top scored docking pose and the ligand bound crystal structure was 0.19 A° (Figure 6). The docked PfDHODH-DMS1 inhibitor complex displayed the same residue-ligand interactions observed in the bound crystal structure demonstrating that the Glide SP docking protocol was capable of regenerating the experimentally observed binding mode of the PfDHODH inhibitor, DSM1.
The initial docking screen yielded 2150 hit compounds; this was narrowed down to 919 compounds by applying the reasonable Glide docking score threshold of ≤ −5 kcal/mol. The compounds were then re-docked into the PfDHODH using the more rigorous docking method, Glide XP, with the same internal docking parameters. Glide XP docking scores for the top pose of each hit compound were in the range of −11.60 to −1.63 kcal/mol. Recent studies have shown that binding energy calculations performed by Prime/MM-GBSA on docked complexes generated by Glide XP can be useful for predicting the binding affinity of hDHODH and PfDHODH inhibitors.35, 63 Prime/MM-GBSA calculations were performed for Glide XP poses with all these poses and based on this, 809 hit compounds, with predicted binding energies ranging from −107.10 to −20.04 kcal/mol, were selected for further analysis.
The predicted binding mode for each of the 809 hit compounds was visually inspected and the type and number of potential interactions within the PfDHODH binding site were analyzed. For this, a set of criteria were derived based on the analyses of critical interactions observed for known malaria specific inhibitors bound to the PfDHODH. 22, 24-26 The inhibitor binding site of PfDHODH is considerably flexible and thus accommodates different classes of inhibitors. The malaria specific trizolopyrimidine-based inhibitors bind to two different binding sites, the hydrogen-bond pocket and the hydrophobic pocket, within the active site of PfDHODH. The trizolopyrimidine ring interacts with the hydrogen-bond pocket, forming hydrogen-bonds between His185 and the bridging nitrogen N1 and between Arg265 and the pyridine nitrogen N5. The ring is also largely surrounded by hydrophobic residues Val532, Leu172, Leu176, Cys184, and Gly181. The naphthyl (DSM1), anthracenyl (DSM2) and phenyl-trifluoromethyl (DSM74) interact with the hydrophobic pocket formed by residues Ile237, Leu189, Leu197, Met536, Phe227, and Phe188. In addition, there is an edge to face stacking interaction between Phe227 and Phe188, which contributes to the potent binding of the inhibitors. An ordered water molecule, W15, which was included in the structure used for docking, mediates a hydrogen-bond between the inhibitors and hydroxyl group of Tyr528. Based on this analysis, the hit compounds were filtered using the following criteria: (i) The compound should have at least one hydrogen-bond interaction with His185, Arg265 and W15 in the hydrogen-bond pocket; (ii) The compound should have a hydrophobic aromatic group facing towards the hydrophobic pocket. Of the 809 hits, 235 compounds satisfied these criteria.
In order to ensure the identification of new and structurally diverse PfDHODH inhibitors and to probe the structural novelty of selected compounds, molecular similarity assessment and cluster analysis were performed for the 235 selected compounds and known PfDHODH inhibitors (Figure S3). Structurally diverse compounds that had a Tanimoto coefficient (Tc) similarity of ≤ 0.6 compared with known PfDHODH inhibitors were selected for further testing. Based on the binding orientation, structural diversity and availability of each sample from the NCI Chemotherapeutic Agents Repository, 39 compounds were initially shipped and tested for in vitro activity against PfDHODH and a chloroquine sensitive strain of P. falciparum (PfNF54) based on the parasite lactate dehydrogenase assay. The topological similarities (Tc) of the selected compounds compared with known inhibitors varied from 0.59 to 0.09 and the average similarity was 0.23 (Tc of 1.0 indicates identical 2D structures). In addition, the selected compounds were compared to known 100 active compounds (Table S2) varied from the maximum similarities from 0.53 to 0.12 and the minimum similarities from 0.09 to 0.03 (Table S3). These analyses show that the selected candidates had new chemotypes.
The 39 virtual hits were evaluated in vitro using a colorimetric PfDHODH enzyme assay that measures DHODH activity by monitoring the reduction of 2, 6-dichloroindophenol. IC50 values were determined for compounds that showed at least 30% PfDHODH inhibition at a concentration of 10 μM in initial inhibition assays. DSM265, a clinical candidate with an IC50 value of 0.008 μM, was used as reference compound. Seven of 39 compounds resulted in 25-65% PfDHODH inhibition at a concentration of 10 μM (Table 4). Compound NSC332161 exhibited 65% inhibition at a concentration of 10 μM and had an IC50 of 5.1 μM. Compound NSC85749 and NSC72405 showed 46 % and 30 % inhibition at a concentration of 10 μM, respectively, and had IC50 of >100 μM. NSC85749 was moderately active at concentrations above 30 μM. The remaining compounds were poor PfDHODH inhibitors displaying <25 % of inhibition at an inhibitor concentration of 10 μM (Table 5 and Tables S4 and S5).
To increase the structural diversity of the compounds, an additional 23 compounds with 2D similarity to NSC332161, NSC85749 and NSC72405 were obtained from the 235 hits of structure-based virtual screening of NCI. The most active inhibitors within the new set of compounds were NSC336047 and NSC343533 which displayed 93% and 50% enzyme inhibition at a concentration of 10 μM and exhibited IC50 values of 0.38 and 20 μM, respectively (Table 4). NSC336047 was evaluated for species-selectivity against hDHODH: it displayed an IC50 value of >100 μM for hDHODH, indicating that it is highly species selective. A771726, a hDHODH specific inhibitor with an IC50 value of 0.45 μM, was used as the reference compound.
Antiplasmodial activities of all selected compounds were evaluated against a Chloroquine sensitive strain of P. falciparum NF54 (PfNF54) based on the parasite lactate dehydrogenase assay (Supporting information) and [3H] hypoxanthine incorporation assay (Supporting information). Chloroquine and Artesunate, potent antimalarial drugs, was used as positive controls in both assays. Chloroquine and Artesunate showed IC50 values of 0.013 μM and 0.010 μM, respectively in the lactate dehydrogenase assay while IC50 values for Chloroquine and Artesunate were 0.012 μM and 0.009 μM, respectively in [3H] hypoxanthine incorporation assay. NSC336047, NSC332161 and NSC343533 showed antiplasmodial activity with IC50 values of 26.28 μM, 23.37 μM and >35.29 μM, respectively (Table 4). In addition, thirteen compounds [nine compounds (Table S4) and four compounds (Table 5)] showed IC50 values of <50 μM against PfNF54 and other compounds showed IC50 values of >50 μM (Table S5). Amongst them, four compounds, which displayed very low activity in PfDHODH inhibition assays, exhibited inhibitory activities against parasite growth in the range of 5-12 μM (Table 5) suggesting that they have a different mechanism of action. Further experimental work is required to elucidate the targets of these 4 compounds.
In order to gain insight into the binding mechanism of the three compounds identified, their binding modes in the active site of PfDHODH were investigated. The predicted binding modes of the three PfDHODH inhibitors are shown in Figure 7. The docking results showed that all three hit compounds adopt very similar orientations within the active site and engage in similar interactions. The dimethyl aniline of NSC336047 (IC50=0.38 μM), the benzodioxole of NSC332161 (IC50=5.1 μM) and the indole of NSC343533 (IC50=20 μM) fit well into the hydrophobic pocket, making favorable hydrophobic interactions with residues Phe188, Leu189, Leu197, Phe227, Ile237, Leu240 and Met536 (Figure 7). In addition, there is likely an edge-to-face π-π stacking interaction between the Phe188/Phe227 and the phenyl ring of these three compounds, which could enhance the potency of the compound.
The hydrophilic part of carboxylic ethyl ester moieties occupies the hydrophilic site of PfDHODH, forming a hydrogen-bond interaction with residue Arg265. In addition, the ester moiety of both NSC336047 and NSC332161 forms hydrogen-bond with a crucial water molecule (W15) (Figure 7). NSC343533 does not interact with W15 as one of its carboxylic ethyl ester moieties is replaced with CN [Figure 7 (E) and (F)]. The bridging nitrogen atom in all three compounds forms an additional hydrogen-bond with His185. Interactions with the His185 have been shown to contribute to the species-selectivity of inhibitors for PfDHODH over hDHODH 26. Additionally, many residues around the inhibitor binding site make hydrophobic contacts with the inhibitors as shown in Figure 7 (B), (D) and (F). The weaker PfDHODH activity observed for NSC343533, compared to the other two active compounds, may be due to the absence of a hydrogen-bond interaction with W15 and the lengthy linker between indole and NH. A comparison of three hit compounds with inactive compounds indicates that the bridging nitrogen is essential for both selectivity and activity for PfDHODH. The inactive molecules may be inactive because they do not adopt the same binding mode as the active compounds.
To explore the structural basis for the observed species selectivity of NSC336047, the compound was docked into the active site of hDHODH (PDB: 1D3G)51 using the Glide XP and Prime/MMGBSA protocols and aligned to the docked PfDHODH-NSC336047 inhibitor complex [Figure 8 (A)]. The binding affinity of NSC336047 to PfDHODH was high with a Glide XP docking score of −9.1 kcal/mol and Prime/MM-GBSA binding energy of −88 kcal/mol (Table 4), while it showed slightly decreased affinity towards hDHODH with a Glide XP docking score of −8 kcal/mol and Prime/MM-GBSA binding energy of −75 kcal/mol, suggesting that the NSC336047 favorably binds to PfDHODH than to hDHODH. The better docking score and binding energy for the PfDHODH-NSC336047 is due to the amino acid replacements in the binding site. As shown in Figure 8 (A), the carboxylic ethyl ester moiety of NSC336047 binds to both PfDHODH and hDHODH in a similar manner. However the dimethyl aniline binds in a completely different position in hDHODH from PfDHODH. The differential binding conformation of dimethyl aniline in hDHODH is because Leu172, Phe188, Leu240 and Met536 in PfDHODH are replaced by Met43, Ala59, Met111 and Pro364 in hDHODH. The bridging nitrogen atom of the compound forms a hydrogen-bond with Ala55 in hDHODH instead of a hydrogen-bond with His185 as observed in the PfDHODH. An additional hydrogen-bond was also observed between Arg136 (corresponding Arg265 in PfDHODH) and the carbonyl group of NSC336047. The comparison suggests that the amino acid replacements, the absence of a hydrogen-bond of His185 and π-stacking interaction of Phe188/Phe227 with NSC336047, provides a basis for the weak inhibitory activity against hDHODH..
In addition, PfDHODH-NSC336047 complex was compared to the crystal structures of PfDHODH bound to the human-specific inhibitor A771726 (PDB ID: 1TV5)50 and malaria-specific inhibitor triazolopyrimidine DSM1 (PDB ID: 3I65)26 [Figure 8 (B)]. The superimposition shows that all three classes of inhibitors nicely overlap in the hydrogen-bond pocket that forms hydrogen-bond interactions with His185 and Arg265. NSC336047 interacts with the hydrophobic pocket in the same manner as the triazolopyrimidine (DSM1) and also forms a likely stacking interaction between NSC336047 and Phe188/Phe227 [Figure 8 (B)], as in the case of DSM1, which has been suggested to contribute to the strong binding of inhibitors.26 However, A771726 occupies a different hydrophobic pocket owing to the altered conformation of Phe188 [Figure 8 (B)], explaining the poor activity of A77 1726 against hDHODH.50
To further investigate the species selectivity of NSC336047, the docked complex of PfDHODH-NSC336047 was aligned to the crystal structure of hDHODH bound to brequinar (PDB ID: 1D3G)51, a hDHODH specific-inhibitor [Figure 8(C)]. The alignment suggests that the selective binding of NSC336047 to PfDHODH is due to the previously mentioned amino acid replacements Met43, Ala59, Met111 and Pro364 in hDHODH for Leu172, Phe188, Leu240 and Met536 in PfDHODH as shown in Figure 8 (C). The replaced residues, Phe188 and Met536 in PfDHODH, block access of brequinar into the binding pocket, thus hindering binding to the PfDHODH structure17, 26 and providing the structural basis for the poor inhibitory activity of brequinar against the malarial parasite enzyme.
The structural similarity of the nine compounds in Table 4, the four compounds in Table 5 and nine compounds (Table S4) identified in this study were analyzed using SciFinder (SciFinder; American Chemical Society, 2015) and compared to known PfDHODH inhibitors in ChEMBL. 53 These databases contain compounds that are annotated based on patents and primary scientific literature. SciFinder indicated that all these compounds have not been reported as either PfDHODH inhibitors or antimalarial compounds. The results of the analysis with ChEMBL showed varying degrees of similarities with Tc values ranging from 0.67 to 0.15 (Supplemental Table. S6). The three PfDHODH inhibitors, NSC336047, NSC332161 and NSC343533, had maximum Tc similarities of 0.67, 0.59 and 0.47, and minimum similarities of 0.04, 0.06 and 0.04, respectively when compared with previously reported PfDHODH inhibitors in ChEMBL (Table S6). While the novelty analysis indicates that these compounds share structural similarity with known PfDHODH inhibitors, these compounds have not previously been reported as PfDHODH inhibitors. The four PfNF54 inhibitors (Table 4) and nine inhibitors (Table S4) had maximum Tc similarities of 0.61 and 0.59, with minimum similarities of 0.03 and 0.03, respectively when compared with ChEMBL PfDHODH compounds (Supplemental Table S6). This analysis shows that these compounds have also not been previously reported as antimalarial compounds.
Over the past decades, there have been several advances in the treatment and control of malaria. Nevertheless, there is still an urgent need for the discovery and development of potent antimalarial drugs that are able to circumvent parasite resistance. The PfDHODH has emerged as a promising selective target for the development of novel therapeutics for malaria as it plays an essential role in de novo pyrimidine biosynthesis. However, limited effort has gone into the development of more effective and selective inhibitors for PfDHODH. To address this we have identified new PfDHODH inhibitors using a combination of pharmacophore and structure-based approaches to screen the NCI database. 3D QSAR pharmacophore models were developed for PfDHODH dihydrothiophenone inhibitors and the top validated model was used as a virtual screening filter. The hit compounds were further prioritized using docking studies and MMGBSA binding free energy calculations. The virtual screening resulted in the identification of nine hit compounds that inhibited PfDHODH by >25% at a concentration of 10 μM, three of which had IC50 ≤ 20 μM. The most active compound showed the species-selective inhibition for PfDHODH (IC50 = 0.38 μM) over hDHODH (IC50 ≥ 100 μM). Thirteen compounds inhibited growth of the Pf NF54 parasite with IC50 values of <50 μM, of which four had IC50 values between 5-12 μM. To the best of our knowledge, these compounds have not been reported previously as PfDHODH and P. falciparum parasite growth inhibitors.
This study demonstrates that the adopted virtual screening strategy is capable of scaffold hopping of already known chemotypes into previously unexplored regions of chemical space. Therefore, the presented virtual screening approach may be used in the hit identification rather than lead optimization. In addition, it could be expected that the systematic SAR analysis of identified hit compounds in the study may be useful for the identification and development of more potent PfDHODH and parasite growth inhibitors.
The authors thank the Drug Synthesis and Chemistry Branch, Developmental Therapeutics Program, Division of Cancer Treatment and Diagnosis of the National Cancer Institute, for kindly providing the compounds screened in this study. We also thank the Centre of High Performance Computing, Cape Town (http://www.chpc.ac.za) for Discovery Studio 4.0 license. This work was supported by funds from the United States National Institutes of Health grant, R01AI103947 (to MAP) and from the Welch Foundation (I-1257) (to MAP). MAP holds the Beatrice and Miguel Elias Distinguished Chair in Biomedical Science and the Carolyn R. Bacon Professorship in Medical Science and Education. The University of Cape Town, South African Medical Research Council, and South African Research Chairs Initiative of the Department of Science and Technology, administered through the South African National Research Foundation are gratefully acknowledged for support (K.C). SW thanks Christoph Fischli for assistance in performing antimalarial assays. EP would like to acknowledge Lauren Arendse for proofreading the manuscript.
Experimental procedure for antiplasmodial activity assays; 2D structures of 19 training set molecules (Figure S1); 2D structures of 19 test set molecules (Figure S2); 2D structures of PfDHODH known inhibitors (Figure S3); Superimposition of ten hypotheses (Hypo1-Hypo10) (Figure S4); Correlation plot for training set and test set compounds (Figure S5); Experimental and predicted activities of the test set compounds by Hypo 1 (Table S1); Active PfDHODH inhibitors used for enrichment factor calculation (Table S2); Comparison of 39 candidates to known 100 active PfDHODH inhibitors (Table S3); Nine hit compounds showing ≤50 μM against PfNF54 (Table S4); Inactive compounds (Table S5); Tanimoto Similarity analysis of hit compounds with ChEMBL database (Table S6). This material is available free of charge via the internet at http://pubs.acs.org.