Protein-Ligand interaction studies on 2, 4, 6- trisubstituted triazine derivatives as anti-malarial DHFR agents using AutoDock.

The dihydrofolate reductase (DHFR) domain of P. falciparum is one of the few well defined targets in malarial chemotherapy. The enzyme catalyzes the nicotinamide adenine dinucleotide phosphate (NADPH) dependent reduction of dihydrofolate to tetrahydrofolate. Protein-ligand interactions were studied using DHFR protein 2BL9, extracted from PDB to evaluate the strength of affinity of various molecules towards ligand binding site and to study the extent of correlation between experimental values and computational dock scores. AutoDock runs resulted in binding energy scores from -7.14 to -10.72 kcal/mol. Among the five inhibitors (Bioorganic and Medicinal Chemistry Letters 15 2005 531-533) selected for docking studies, an excellent correlation was observed in all cases, for instance, experimentally reported most active molecule 2a (MIC: 1µg/ml) showed a high dock score (-10.72 kcal/mol) than the remaining inhibitors. Therefore, molecular docking using AutoDock suggests the importance of evaluating the prediction accuracy of various molecules as evidenced by a correlation coefficient of 0.961 between experimental activities and AutoDock binding energies.

Background

Dihydrofolate Reductase is a small enzyme that plays a supporting, but an essential role, in the building of DNA and other processes. DHFR catalyzes the reduction of 7, 8-dihydrofolate (DHF) through stereo-specific hydride transfer from reduced nicotinamide-adenine dinucleotide phosphate (NADPH) cofactor. The enzyme is essential for tetrahydrofolate (THF) biosynthesis, plays a central role in promoting cell growth and proliferation, and is the target of several anticancer and antibiotic drugs [1]. Enzymes with essential roles are sensitive targets for drug therapy. Dihydrofolate reductase was the first enzyme to be targeted for cancer chemotherapy [2]. PfDHFR is a well-characterized target for antimalarial antifolate drugs, such as pyrimethamine (Pyr) and cycloguanil (Cyc). DHFR of Plasmodium and other protozoa are expressed on the same polypeptide chain with its accompanying enzyme, thymidylate synthase (TS). TS and DHFR sequentially catalyze consecutive reactions in the biosynthesis of thymidylate, which is required for DNA synthesis [3].

Malaria, caused by the protozoan parasite Plasmodium falciparum, affects about 300 million individuals and causes about 2 million deaths per year. Traditional anti malarial agents such as chloroquine are ineffective in many regions of the world due to drug resistance, chemotherapy targeted at dihydrofolate Reductase (DHFR) and thymidylate synthase (TS) has proven to be highly effective. Even partial inhibition of DHFR or TS can lead to DNA strand fragmentation and cell death [4– 5]. In malaria pharmacology, these two enzymes are of particular interest because traditional drugs such as pyrimethamine and cycloguanil are known to inhibit parasite DHFR-TS. In recent years, however, the effectiveness of these inhibitors has been compromised by malarial parasite strains expressing mutant forms of DHFRTS [6]. One obvious way to continue selective killing of malarial parasites is to identify new folate analogs directed at DHFR- or TS-active sites that are effective against drug-resistant parasites [7]. The folate antagonists are an important class of therapeutic compounds, as evidenced by their use as antiinfective, anti neoplastic, and anti-inflammatory drugs [8]. Dihydrofolate reductase -thymidylate synthase (DHFR-TS) from Plasmodium falciparum is a validated target for antifolate antimalarials [9]. Hence, inhibition of DHFR and its role in combating malaria has provided us the rationale to carry out structure based drug design studies. In this study, compounds synthesized to target malaria, from literature, are considered to perform DHFR protein, extracted from Protein Data Bank, and ligand interactions using docking software, AutoDock. The validation procedures are carried out and the interactions and scores are generated.

Materials and Methodology

Out of all the entries for DHFR (dihydrofolate reductase) from RCSB protein data bank, Plasmodium vivax 2BL9 was taken for docking analysis. The active site residues were found to be THR194, ASP53, MET54, PHE57, LEU45, SER117, ALA15, ILE13 and ILE173. Protein - ligand docking studies were performed to evaluate the algorithm and scoring function efficiency between a standalone AutoDock 3.0.5 and experimental activities. The five 2, 4, 6- trisubstituted triazine inhibitors (active, moderately active and inactive) were selected from the article Bioorganic and Medicinal Chemistry Letters 15 2005 531-533 for docking studies. All these molecules as well as the bound ligand of the protein 2BL9 were docked by using the software AutoDock and the binding affinities were predicted. All molecules are drawn using ISIS draw and energy minimized using Tsar Software (www.accelrys.com).

Automated docking was used to locate the appropriate binding orientations and conformations of various inhibitors in the 2BL9 binding pocket. To perform the task, genetic algorithm routine implemented in the program AutoDock 3.0. was employed [10]. All water molecules were removed from the original Protein Data Bank file. Polar hydrogen atoms were added and Kolllman charge, atomic solvation parameters and fragmental volumes were assigned to the protein using AutoDock Tools (ADT). To validate the docking protocol, bound ligand CP61240 coordinates in the crystal complex was removed and the bond orders were checked. For docking calculations, Gasteiger partial charges were assigned to the tested derivatives and CP61240, and non-polar hydrogen atoms were merged. All torsions were allowed to rotate during docking.

The program AutoGrid used to generate the grid maps. Each grid was centered at the crystal structure of the corresponding 2BL9 bound ligand CP61240. The grid dimensions were 60 × 60 × 60 A°3 with points separated by 0.375 A°. Lennard-Jones parameters 12-10 and 12-6, supplied with the program, were used for modeling H-bonds and van der Waals interactions, respectively. The distance-dependent dielectric permittivity of Mehler and Solmajer [11] was used for calculation of the electrostatic grid maps. For all ligands, random starting positions, random orientations and torsions were used. The translation, quaternion and torsion steps were taken from default values in AutoDock. The Lamarckian genetic algorithm and the pseudo-Solis and Wets methods were applied for minimization using default parameters. The standard docking protocol for rigid and flexible ligand docking consisted of 10 independent runs per ligand, using an initial population of 50 randomly placed individuals, with 2.5 × 106 energy evaluations, a maximum number of 27000 iterations, a mutation rate of 0.02, a crossover rate of 0.80, and an elitism value of 1. The probability of performing a local search on an individual in the population was 0.06, using a maximum of 300 iterations per local search. After docking, the 10 solutions were clustered into groups with RMS deviations lower than 1.0 Å. The clusters were ranked by the lowest energy representative of each cluster.

Results and Discussion

Before docking the ligands (Table 1 see Table 1) into the 2BL9 binding site, the docking protocol was validated. CP61240 bound ligand was removed from the active site and docked back into the binding pocket. AutoDock predicted binding conformation of CP61240 is shown in Figure 1 with the X-ray crystallographic obtained conformational superposition [12]. The RMSD of all atoms between these two conformations is 0.72 A° indicating that the parameters for the docking simulation are reasonable in reproducing the X-ray crystal structure and can be extended to search the enzyme binding conformations for other inhibitors accordingly. The ligands selected for docking studies and the corresponding interaction energies are given in Table 1. The active site of 2BL9 offers many different binding modes for these compounds as they are strongly dependent on the attached substituent. This observation results from the comparison of ligand binding orientations in the active site for known CP61240 and selected compounds. All the ligands were docked deeply within the binding pocket region of 2BL9 forming more hydrogen bonds with Ile173, Asp53 and Cys14 (Figure 2 and Figure 3). As shown in Table 1, their AutoDock binding free energies (ΔGb, kcal/mol) and inhibition constants were obtained. Among them compound1 exhibited the lowest free energy -7.14 kcal/mol. In other words, it posses the highest potential binding affinity into the binding site of 2BL9 protein. The overall good correlation, given in Figure 4, between the experimental biological activities and the AutoDock binding affinities made a clear representation of utilizing computational techniques for drug design. The correlation coefficient of 0.961 indicates that the correlation was fairly good and demonstrates the applicability of the method.

Figure 1

Docked superimposed image of 2BL9 bound ligand CP61240 within the binding site region.

Figure 2

Image showing compound 2a interacting residues within the active site region of 2BL9.

Figure 3

Image showing Compound 1 interacting residues within the active site region of 2Bl9.

Figure 4

Correlation between experimental and AutoDock binding free energies of selected 2, 4, 6-trisubstituted triazine ligands.

Most docked inhibitors interacted in the same fashion as observed with cocrystallized ligand of 2BL9 protein showing H-bonds with Ile13, Asp53, Thr194 and Ile173. The binding modes and geometrical orientation of compounds 2a, 2b, 2c, 2l and compound 1 nearly resemble the 2BL9 bound ligand CP61240. Moreover, the compound 2a showed two hydrogen bonds with Ala15, and the compound 2b showed only one interaction with Cys49 respectively, where as compound 2c has no H-bond interactions. Compounds 2l and 1, showed low binding energy than others, yet the number of hydrogen bond interactions are two with Thr194 and Ser117 and four with Thr194,Ile13,Asp53 and Ile173, respectively. The most active molecule was found to be 2a and the inactive compound was compound 1. Compound 2a exhibited the dock score -10.72 kcal/mol and it showed two H-bond interactions (Figure 2) with active site of 2BL9, one with carbonyl oxygen of Ala15 and the other with nitrogen of Ala15. From the Table 1, it became evident that there existed a good correlation between experimental activities and computational binding energies as obtained with AutoDock. Hence, further proof was provided by plotting a graph between experimental values and dock scores in Figure 4, where it is clear that they represented a correlation of 0.961.

Conclusion

AutoDock runs resulted in binding energy scores that range from -7.14 to- 10.72 kcal/mol. Among the five 2, 4, 6-trisubstituted triazine inhibitors selected for docking studies, an excellent correlation of 0.96 was observed. Therefore, the dock analysis suggests the importance of evaluating the prediction accuracy of scoring functions adopted in AutoDock and showed excellent correlation between binding free energy (kcal/mol) and IC50 (µM) values against 2BL9 DHFR protein. Hence, this study demonstrates the importance of drug design studies and utility of computational tools in correlating experimental values with computational binding energy scores.