Novel approach to structure-based pharmacophore search using computational geometry and shape matching techniques

Computationally efficient structure-based virtual screening methods have recently been reported that seek to find effective means to utilize experimental structure information without employing detailed molecular docking calculations. These tools can be coupled with efficient experimental screening technologies to improve the probability of identifying hits and leads for drug discovery research. Commercial software ROCS (rapid overlay of chemical structures) from Open Eye Scientific is such an example, which is a shape-based virtual screening method using the 3D structure of a ligand, typically from a bound X-ray costructure, as the query. We report here the development of a new structure-based pharmacophore search method (called Shape4) for virtual screening. This method adopts a variant of the ROCS shape technology and expands its use to work with an empty crystal structure. It employs a rigorous computational geometry method and a deterministic geometric casting algorithm to derive the negative image (i.e., pseudoligand) of a target binding site. Once the negative image (or pseudoligand) is generated, an efficient shape comparison algorithm in the commercial OE SHAPE Toolkit is adopted to compare and match small organic molecules with the shape of the pseudoligand. We report the detailed computational protocol and its computational validation using known biologically active compounds extracted from the WOMBAT database. Models derived for five selected targets were used to perform the virtual screening experiments to obtain the enrichment data for various virtual screening methods. It was found that our approach afforded similar or better enrichment ratios than other related methods, often with better diversity among the top ranking computational hits.     

Using consensus-shape clustering to identify promiscuous ligands and protein targets and to choose the right query for shape-based virtual screening

Ligand-based shape matching approaches have become established as important and popular virtual screening (VS) techniques. However, despite their relative success, many authors have discussed how best to choose the initial query compounds and which of their conformations should be used. Furthermore, it is increasingly the case that pharmaceutical companies have multiple ligands for a given target and these may bind in different ways to the same pocket. Conversely, a given ligand can sometimes bind to multiple targets, and this is clearly of great importance when considering drug side-effects. We recently introduced the notion of spherical harmonic-based "consensus shapes" to help deal with these questions. Here, we apply a consensus shape clustering approach to the 40 protein-ligand targets in the DUD data set using PARASURF/PARAFIT. Results from clustering show that in some cases the ligands for a given target are split into two subgroups which could suggest they bind to different subsites of the same target. In other cases, our clustering approach sometimes groups together ligands from different targets, and this suggests that those ligands could bind to the same targets. Hence spherical harmonic-based clustering can rapidly give cross-docking information while avoiding the expense of performing all-against-all docking calculations. We also report on the effect of the query conformation on the performance of shape-based screening of the DUD data set and the potential gain in screening performance by using consensus shapes calculated in different ways. We provide details of our analysis of shape-based screening using both PARASURF/PARAFIT and ROCS, and we compare the results obtained with shape-based and conventional docking approaches using MSSH/SHEF and GOLD. The utility of each type of query is analyzed using commonly reported statistics such as enrichment factors (EF) and receiver-operator-characteristic (ROC) plots as well as other early performance metrics.     

Combined vibrational,conformational energy and electron diffraction studies of trans-2-decalone

The structure of trans-2-decalone was investigated by combined gas phase electron diffraction, conformational energy and vibrational analyses. In this study first the minimum energy conformations for trans-2-decalone were calculated by molecular mechanics techniques using a force field described in the literature; the same field and the minimum energy conformations were then used in subsequent vibrational analyses to calculate the mean amplitudes of vibration for each minimum energy conformation of trans-2-decalone; these mean amplitudes and the corresponding internuclear distances were then used to calculate theoretical electron diffraction radial distribution curves which were compared to the experimental curves. Four conformational energy minima were investigated with one or both rings in a chair or non-chair form. The results of the combined investigations indicate that the molecule exists in the conformation which has both rings in a distorted chair form.     

ROCS-derived features for virtual screening

Rapid overlay of chemical structures (ROCS) is a standard tool for the calculation of 3D shape and chemical (“color”) similarity. ROCS uses unweighted sums to combine many aspects of similarity, yielding parameter-free models for virtual screening. In this report, we decompose the ROCS color force field into color components and color atom overlaps, novel color similarity features that can be weighted in a system-specific manner by machine learning algorithms. In cross-validation experiments, these additional features significantly improve virtual screening performance relative to standard ROCS.     

The MPSim-Dock hierarchical docking algorithm: application to the eight trypsin inhibitor cocrystals

To help improve the accuracy of protein-ligand docking as a useful tool for drug discovery, we developed MPSim-Dock, which ensures a comprehensive sampling of diverse families of ligand conformations in the binding region followed by an enrichment of the good energy scoring families so that the energy scores of the sampled conformations can be reliably used to select the best conformation of the ligand. This combines elements of DOCK4.0 with molecular dynamics (MD) methods available in the software, MPSim. We test here the efficacy of MPSim-Dock to predict the 64 protein-ligand combinations formed by starting with eight trypsin cocrystals, and crossdocking the other seven ligands to each protein conformation. We consider this as a model for how well the method would work for one given target protein structure. Using as a criterion that the structures within 2 kcal/mol of the top scoring include a conformation within a coordinate root mean square (CRMS) of 1 A of the crystal structure, we find that 100% of the 64 cases are predicted correctly. This indicates that MPSim-Dock can be used reliably to identify strongly binding ligands, making it useful for virtual ligand screening.     

A robust force field based method for calculating conformational energies of charged drug-like molecules

The binding affinity of a drug-like molecule depends among other things on the availability of the bioactive conformation. If the bioactive conformation has a significantly higher energy than the global minimum energy conformation, then the molecule is unlikely to bind to its target. Determination of the global minimum energy conformation and calculation of conformational penalties of binding is a prerequisite for prediction of reliable binding affinities. Here, we present a simple and computationally efficient procedure to estimate the global energy minimum for a wide variety of structurally diverse molecules, including polar and charged compounds. Identifying global energy minimum conformations of such compounds with force field methods is problematic due to the exaggeration of intramolecular electrostatic interactions. We demonstrate that the global energy minimum conformations of zwitterionic compounds generated by conformational analysis with modified electrostatics are good approximations of the conformational distributions predicted by experimental data and with molecular dynamics performed in explicit solvent. Finally the method is used to calculate conformational penalties for zwitterionic GluA2 agonists and to filter false positives from a docking study.     

Exploring polypharmacology using a ROCS-based target fishing approach

Polypharmacology has emerged as a new theme in drug discovery. In this paper, we studied polypharmacology using a ligand-based target fishing (LBTF) protocol. To implement the protocol, we first generated a chemogenomic database that links individual protein targets with a specified set of drugs or target representatives. Target profiles were then generated for a given query molecule by computing maximal shape/chemistry overlap between the query molecule and the drug sets assigned to each protein target. The overlap was computed using the program ROCS (Rapid Overlay of Chemical Structures). We validated this approach using the Directory of Useful Decoys (DUD). DUD contains 2950 active compounds, each with 36 property-matched decoys, against 40 protein targets. We chose a set of known drugs to represent each DUD target, and we carried out ligand-based virtual screens using data sets of DUD actives seeded into DUD decoys for each target. We computed Receiver Operator Characteristic (ROC) curves and associated area under the curve (AUC) values. For the majority of targets studied, the AUC values were significantly better than for the case of a random selection of compounds. In a second test, the method successfully identified off-targets for drugs such as rimantadine, propranolol, and domperidone that were consistent with those identified by recent experiments. The results from our ROCS-based target fishing approach are promising and have potential application in drug repurposing for single and multiple targets, identifying targets for orphan compounds, and adverse effect prediction.     

Do the molecules which form discotic liquid crystals have disc-like structures? The conformation of a simple model compound, 1,2-dihydroxydiacetylbenzene, determined from the NMR spectra of samples dissolved in liquid crystalline solvents

Many discotic mesogens are molecules with a central aromatic ring with adjacent alkylcarboxylate substituents. The simplest such molecule, 1,2-dihydroxydiacetylbenzene, which is not mesogenic, is studied by NMR spectroscopy as a solute in a nematic solvent. The spectra are analysed to give sets of residual dipolar couplings, D_ij , which are then used to test models for the conformation adopted by the acetate side groups. The conformations and geometry of an isolated molecule are calculated by the ab initio MP2/6-311G method and also by the DFT approach using the B3LYP functional with the 6-311++G** basis set. The quantum chemical calculations find that the minimum energy conformer has the acetate groups rotated in opposite directions out of the ring plane, and this kind of structure is also consistent with the NMR data.     

The histamine H1-receptor antagonist binding site. Part I: Active conformation of cyproheptadine

Summary The active conformation of several histamine H₁-antagonists is investigated. As a template molecule we used the antagonist cyproheptadine, which consists of a piperidylene ring connected to a tricyclic system. The piperidylene moiety is shown to be flexible. The global minimum is a chair conformation but, additionally, a second chair and various boat conformations have to be considered, as their energies are less than 5 kcal/mol above the energy of the global minimum. Two semi-rigid histamine H₁-antagonists, phenindamine and triprolidine, were fitted onto the various conformations of cyproheptadine in order to derive the pharmacologically active conformation of cyproheptadine. At the same time, the active conformation of both phenindamine and triprolidine was derived. It is demonstrated that, within the receptor-bound conformation of cyproheptadine, the piperidylene ring most probably exists in a boat form.     

Kinetics of cooperative protein folding involving two separate conformational families

The master equation that describes the kinetics of protein folding is solved numerically for a portion of Staphylococcal Protein A by a Laplace transformation. The calculations are carried out with 50 local-minimum conformations belonging to two conformational families. The master equation allows for transitions among all the 50 conformations in the evolution toward the final folded equilibrium distribution of conformations. It is concluded that the native protein folds in a fast cooperative process. The global energy minimum of a native protein can be reached after a sufficiently long folding time regardless of the initial state and the existence of a large number of local energy minima. Conformations representing non-native states of the protein can transform to the native state even if they do not belong to the native conformational family. Given a starting conformation, the protein molecule can fold to its final conformation through different paths. Finally, when the folding reaches the equilibrium distribution, the protein molecule adopts a set of conformations in which the global minimum has the largest average probability.     

Representing receptor flexibility in ligand docking through relevant normal modes

Inspired by the current representation of the ligand-receptor binding process, a normal-mode-based methodology is presented to incorporate receptor flexibility in ligand docking and virtual screening. However, the systematic representation of the deformation space grows geometrically with the number of modes, and furthermore, midscale loop rearrangements like those found in protein kinase binding pockets cannot be accounted for with the first lowest-frequency modes. We thus introduced a measure of relevance of normal modes on a given region of interest and showed that only very few modes in the low-frequency range are necessary and sufficient to describe loop flexibility in cAMP-dependent protein kinase. We used this approach to generate an ensemble of representative receptor backbone conformations by perturbing the structure along a combination of relevant modes. Each ensemble conformation is complexed with known non-native binders to optimize the position of the binding-pocket side chains through a full flexible docking procedure. The multiple receptor conformations thus obtained are used in a small-scale virtual screening using receptor ensemble docking. We evaluated this algorithm on holo and apo structures of cAMP-dependent protein kinase that exhibit backbone rearrangements on two independent loop regions close to the binding pocket. Docking accuracy is improved, since the ligands considered in the virtual screening docked within 1.5 A to at least one of the structures. The discrimination between binders and nonbinders is also enhanced, as shown by the improvement of the enrichment factor. This constitutes a new step toward the systematic integration of flexible ligand-flexible receptor docking tools in structure-based drug discovery.     

Docking and scoring with ICM: the benchmarking results and strategies for improvement

Flexible docking and scoring using the internal coordinate mechanics software (ICM) was benchmarked for ligand binding mode prediction against the 85 co-crystal structures in the modified Astex data set. The ICM virtual ligand screening was tested against the 40 DUD target benchmarks and 11-target WOMBAT sets. The self-docking accuracy was evaluated for the top 1 and top 3 scoring poses at each ligand binding site with near native conformations below 2?? RMSD found in 91 and 95% of the predictions, respectively. The virtual ligand screening using single rigid pocket conformations provided the median area under the ROC curves equal to 69.4 with 22.0% true positives recovered at 2% false positive rate. Significant improvements up to ROC AUC?=?82.2 and ROC((2%))?=?45.2 were achieved following our best practices for flexible pocket refinement and out-of-pocket binding rescore. The virtual screening can be further improved by considering multiple conformations of the target.     

Novel inhibitor discovery through virtual screening against multiple protein conformations generated via ligand-directed modeling: a maternal embryonic leucine zipper kinase example

Kinase targets have been demonstrated to undergo major conformational reorganization upon ligand binding. Such protein conformational plasticity remains a significant challenge in structure-based virtual screening methodology and may be approximated by screening against an ensemble of diverse protein conformations. Maternal embryonic leucine zipper kinase (MELK), a member of serine-threonine kinase family, has been recently found to be involved in the tumerogenic state of glioblastoma, breast, ovarian, and colon cancers. We therefore modeled several conformers of MELK utilizing the available chemogenomic and crystallographic data of homologous kinases. We carried out docking pose prediction and virtual screening enrichment studies with these conformers. The performances of the ensembles were evaluated by their ability to reproduce known inhibitor bioactive conformations and to efficiently recover known active compounds early in the virtual screen when seeded with decoy sets. A few of the individual MELK conformers performed satisfactorily in reproducing the native protein-ligand pharmacophoric interactions up to 50% of the cases. By selecting an ensemble of a few representative conformational states, most of the known inhibitor binding poses could be rationalized. For example, a four conformer ensemble is able to recover 95% of the studied actives, especially with imperfect scoring function(s). The virtual screening enrichment varied considerably among different MELK conformers. Enrichment appears to improve by selection of a proper protein conformation. For example, several holo and unliganded active conformations are better to accommodate diverse chemotypes than ATP-bound conformer. These results prove that using an ensemble of diverse conformations could give a better performance. Applying this approach, we were able to screen a commercially available library of half a million compounds against three conformers to discover three novel inhibitors of MELK, one from each template. Among the three compounds validated via experimental enzyme inhibition assays, one is relatively potent (15; K(d) = 0.37 μM), one moderately active (12; K(d) = 3.2 μM), and one weak but very selective (9; K(d) = 18 μM). These novel hits may be utilized to assist in the development of small molecule therapeutic agents useful in diseases caused by deregulated MELK, and perhaps more importantly, the approach demonstrates the advantages of choosing an appropriate ensemble of a few conformers in pursuing compound potency, selectivity, and novel chemotypes over using single target conformation for structure-based drug design in general.     

Powerful simulated-annealing algorithm locates global minimum of protein-folding potentials from multiple starting conformations

Protein-folding potentials, designed with the explicit goal that the global energy minimum correspond to crystallographically observed conformations of protein molecules, may offer great promise toward calculating native protein structures. Achieving this promise, however, depends on finding an effective means of dealing with the multiple-minimum problem inherent in such potentials. In this study, a protein-folding-potential test system has been developed that exhibits the properties of general protein-folding potentials yet has a unique well-defined global energy minimum corresponding to the crystallographically determined conformation of the test molecule. A simulated-annealing algorithm is developed that locates the global minimum of this potential in four of eight test runs from random starting conformations. Exploration of the energy-conformation surface of the potential indicates that it contains the numerous local minima typical of protein-folding potentials and that the global minimum is not easily located by conventional minimization procedures. When the annealing algorithm is applied to a previously developed actual folding potential to analyze the conformation of avian pancreatic polypeptide, a new conformer is located that is lower in energy than any conformer located in previous studies using a variety of minimization techniques.     

Classifying the conformations of a chemical system using matrices of integers

A chemical system is a collection of atomic nuclei and accompanying electrons interacting by electromagnetic forces and organized into one or more molecules, ions, and transient structures. A general problem is to devise a mathematical representation of the conformation of a system which is useful in further mechanical analysis. Through dihedral rotations of bonds and relative motion of structures, many chemical systems can attain an infinite number of conformations. This paper is a theoretical presentation of a scheme to partition these conformations into a finite number of sets. The scheme involves a discrete aspect of chemical kinematics and uses matrices of integers called proximity matrices. The partitioning scheme is anticipated to be useful in studying reaction mechanisms and interactions between molecules, and in finding conformations of particular interest such as that with a potential energy near the global minimum for a molecule.     

Combined electron diffraction,conformational energy and vibrational investigations of cis-1,4-ditertiary butyl-cyclohexane

The structure of cis-1,4-ditertiarybutylcyclohexane(DTBC) was investigated by combined electron diffraction, conformational and vibrational analyses in order to obtain results which are more conclusive than those previously obtained by electron diffraction alone. In this study, first the minimum energy conformations for DTBC were calculated by the Westheimer-Hendrickson procedure using various force fields described in the literature; the same fields and the minimum energy conformations were used in subsequent vibrational analyses to calculate the mean amplitudes of vibration for each minimum energy conformation of DTBC; these mean amplitudes and the corresponding internuclear distances were then used to calculate the theoretical electron diffraction radial distribution curves which were compared to the experimental curves. The results indicate that the conformational energies of all the minimum energy chair and non-chair forms of DTBC are very similar. In excellent agreement with this, the theoretical radial distribution curves of all minimum energy forms have to be mixed for a best fit to the experimental radial distribution curve. A least squares analysis of the mixture under the described conditions yields for 110 °C a composition of approximately one third chair and two thirds non-chair forms. The quality of the empirical conformational force fields has a definite influence on the reliability of these results.     

Molecular recognition in cyclodextrin complexes of amino acid derivatives: the effects of kinetic energy on the molecular recognition of a pseudopeptide in a nonconstraining host environment as revealed by a temperature-dependent crystallographic study

The crystal structure of a triclinic 2:2 inclusion complex of beta-cyclodextrin with N-acetyl-L-phenylalanine methyl ester has been determined at several temperatures between 298 and 20 K to further study molecular recognition using solid-state supramolecular beta-cyclodextrin complexes. The study reveals kinetic energy dependent changes in guest molecule conformations, orientations, and positions in the binding pocket presented by the crystal lattice. Accompanying these changes are observable differences in guest-guest interactions and hydrogen-bonding interactions in the binding pocket that involve guest molecules, water of hydration molecules, and beta-cyclodextrin molecules. On the basis of the differences observed in the crystal structures, we present a solid-state example of a system that displays the properties of both classical and quantum chemical models. At higher temperatures, the structure conforms to a classical mechanical model with dynamic disorder. At lower temperatures, the observations conform to examples in which there is static disorder representative of models in which quantum states differing in conformation, position, and orientation of components in the crystal structure are occupied. Ab initio theoretical calculations on the different guest molecule conformations have been carried out. Superpositions of theoretical electrostatic surface potential diagrams on the observed molecular positions in the complexes provide confidence that the deconvolution of the guest molecule disorder is acceptable. Temperature-dependent solid-state magic angle spinning deuteron NMR measurements provide evidence for large-amplitude, diffusive motion on a microsecond time scale in the complex.