Novel statistical-thermodynamic methods to predict protein-ligand binding positions using probability distribution functions

We present two novel methods to predict native protein-ligand binding positions. Both methods identify the native binding position as the most probable position corresponding to a maximum of a probability distribution function (PDF) of possible binding positions in a protein active site. Possible binding positions are the origins of clusters composed, on the basis of root-mean square deviations (RMSD), from the multiple ligand positions determined by a docking algorithm. The difference between the methods lies in the ways the PDF is derived. To validate the suggested methods, we compare the averaged RMSD of the predicted ligand docked positions relative to the experimentally determined positions for a set of 135 PDB protein-ligand complexes. We demonstrate that the suggested methods improve docking accuracy by as much as 21-24% in comparison with a method that simply identifies the binding position as the energy top-scored ligand position.     

Predicting Ligand Binding Affinity with Alchemical Free Energy Methods in a Polar Model Binding Site

We present a combined experimental and modeling study of organic ligand molecules binding to a slightly polar engineered cavity site in T4 lysozyme (L99A/M102Q). For modeling, we computed alchemical absolute binding free energies. These were blind tests performed prospectively on 13 diverse, previously untested candidate ligand molecules. We predicted that eight compounds would bind to the cavity and five would not; 11 of 13 predictions were correct at this level. The RMS error to the measurable absolute binding energies was 1.8 kcal/mol. In addition, we computed “relative” binding free energies for six phenol derivatives starting from two known ligands: phenol and catechol. The average RMS error in the relative free energy prediction was 2.5 kcal/mol (phenol) and 1.1 kcal/mol (catechol). To understand these results at atomic resolution, we obtained x-ray co-complex structures for nine of the diverse ligands and for all six phenol analogs. The average RMSD of the predicted pose to the experiment was 2.0 Å (diverse set), 1.8 Å (phenol-derived predictions), and 1.2 Å (catechol-derived predictions). We found that predicting accurate affinities and rank-orderings required near-native starting orientations of the ligand in the binding site. Unanticipated binding modes, multiple ligand binding, and protein conformational change all proved challenging for the free energy methods. We believe that these results can help guide future improvements in physics-based absolute binding free energy methods.     

D2N: Distance to the native

Root-mean-square-deviation (RMSD), of computationally-derived protein structures from experimentally determined structures, is a critical index to assessing protein-structure-prediction-algorithms (PSPAs). The development of PSPAs to obtain 0 Å RMSD from native structures is considered central to computational biology. However, till date it has been quite challenging to measure how far a predicted protein structure is from its native — in the absence of a known experimental/native structure. In this work, we report the development of a metric “D2N” (distance to the native) — that predicts the “RMSD” of any structure without actually knowing the native structure. By combining physico-chemical properties and known universalities in spatial organization of soluble proteins to develop D2N, we demonstrate the ability to predict the distance of a proposed structure to within ± 1.5 ? error with a remarkable average accuracy of 93.6% for structures below 5 ? from the native. We believe that this work opens up a completely new avenue towards assigning reliable structures to whole proteomes even in the absence of experimentally determined native structures. The D2N tool is freely available at http://www.scfbio-iitd.res.in/software/d2n.jsp.     

Ligand-supported homology modelling of protein binding-sites using knowledge-based potentials

A new approach, MOBILE, is presented that models protein binding-sites including bound ligand molecules as restraints. Initially generated, homology models of the target protein are refined iteratively by including information about bioactive ligands as spatial restraints and optimising the mutual interactions between the ligands and the binding-sites. Thus optimised models can be used for structure-based drug design and virtual screening. In a first step, ligands are docked into an averaged ensemble of crude homology models of the target protein. In the next step, improved homology models are generated, considering explicitly the previously placed ligands by defining restraints between protein and ligand atoms. These restraints are expressed in terms of knowledge-based distance-dependent pair potentials, which were compiled from crystallographically determined protein-ligand complexes. Subsequently, the most favourable models are selected by ranking the interactions between the ligands and the generated pockets using these potentials. Final models are obtained by selecting the best-ranked side-chain conformers from various models, followed by an energy optimisation of the entire complex using a common force-field. Application of the knowledge-based pair potentials proved efficient to restrain the homology modelling process and to score and optimise the modelled protein-ligand complexes. For a test set of 46 protein-ligand complexes, taken from the Protein Data Bank (PDB), the success rate of producing near-native binding-site geometries (rmsd<2.0A) with MODELLER is 70% when the ligand restrains the homology modelling process in its native orientation. Scoring these complexes with the knowledge-based potentials, in 66% of the cases a pose with rmsd <2.0A is found on rank 1. Finally, MOBILE has been applied to two case studies modelling factor Xa based on trypsin and aldose reductase based on aldehyde reductase.     

Selection and flexible optimization of binding modes from conformation ensembles

This paper describes a new semi-flexible docking approach named Fado (flexible alignment and docking), which incorporates flexibility by using an ensemble of precomputed ligand conformers. A primary ligand is defined as a linear combination over all input conformers. An optimization with regard to the linear coefficients makes the ligand flexible. Initially, a point matching problem utilizing the Merck Molecular Force Field (MMFF) is modeled in order to compute the correct orientation of the ligand with respect to the target. The problem is then solved through a local optimization approach (RPROP). This is done for 20 randomized ligand orientations, yielding 20 binding modes per complex. Evaluating these modes illustrates that our method is able to reproduce the binding modes of molecules within a few minutes of CPU time. A representative dataset of diverse protein-ligand complexes could be reproduced with 78% accuracy below 2A RMSD distance to the reference crystal structure. Fado is available upon request to the authors (see also http://www.zib.eu/Numerik/projects/docking/projectlong.en.html).     

Prediction of ligand-binding sites of proteins by molecular docking calculation for a random ligand library

A new approach to predicting the ligand-binding sites of proteins was developed, using protein-ligand docking computation. In this method, many compounds in a random library are docked onto the whole protein surface. We assumed that the true ligand-binding site would exhibit stronger affinity to the compounds in the random library than the other sites, even if the random library did not include the ligand corresponding to the true binding site. We also assumed that the affinity of the true ligand-binding site would be correlated to the docking scores of the compounds in the random library, if the ligand-binding site was correctly predicted. We call this method the molecular-docking binding-site finding (MolSite) method. The MolSite method was applied to 89 known protein-ligand complex structures extracted from the Protein Data Bank, and it predicted the correct binding sites with about 80-99% accuracy, when only the single top-ranked site was adopted. In addition, the average docking score was weakly correlated to the experimental protein-ligand binding free energy, with a correlation coefficient of 0.44.     

A Method for Helical RNA Global Structure Determination in Solution Using Small-Angle X-Ray Scattering and NMR Measurements

We report a “top-down” method that uses mainly duplexes' global orientations and overall molecular dimension and shape restraints, which were extracted from experimental NMR and small-angle X-ray scattering data, respectively, to determine global architectures of RNA molecules consisting of mostly A-form-like duplexes. The method is implemented in the G2G (from global measurement to global structure) toolkit of programs. We demonstrate the efficiency and accuracy of the method by determining the global structure of a 71-nt RNA using experimental data. The backbone root-mean-square deviation of the ensemble of the calculated global structures relative to the X-ray crystal structure is 3.0 ± 0.3 Å using the experimental data and is only 2.5 ± 0.2 Å for the three duplexes that were orientation restrained during the calculation. The global structure simplifies interpretation of multidimensional nuclear Overhauser spectra for high-resolution structure determination. The potential general application of the method for RNA structure determination is discussed.     

First report on exploring structural requirements of alpha and beta thymidine analogs for PfTMPK inhibitory activity using in silico studies

With the emergence of multi-drug resistance of the currently available antimalarial drugs including the “magic bullet” artemisinin derivatives in the market, there is an urgent need for discovery and development of new potent antimalarial molecules. The present work deals with quantitative structure–activity relationship (QSAR) modeling, pharmacophore mapping and docking studies of a series of 35 thymidine analogs as inhibitors of Plasmodium falciparum thymidylate kinase (PfTMPK), an enzyme that catalyzes phosphorylation of thymidine monophosphate (TMP) to thymidine diphosphate (TDP). The models were validated both internally and externally and significant statistical results were obtained, indicating the robustness and reliability of the developed models. The docking study was performed using the LigandFit option of receptor–ligand interactions protocol section available in Discovery Studio 2.1 where lower RMSD values (0.6931 Å) between the co-crystallized ligand and re-docked ligand assured that the ligand was bound in the same binding pocket. The QSAR, pharmacophore mapping and docking studies provide an understanding of important structural requirements or essential molecular properties, or features of molecules, and important binding interactions, and provide an important guidance for the chemist to synthesis of new molecules with improved PfTMPK inhibitory activity profile. This work revealed the importance of –NH-fragment, electrophilicity of the molecules and the number of oxygen atom towards the PfTMPK inhibitory activity of the molecules. To the best of our knowledge, this work presents the first QSAR and pharmacophore report for thymidine analogs which may serve as an efficient tool for the design and synthesis of potent molecules as PfTMPK inhibitors to address the increasing threat of multi-drug resistance against P. falciparum.     

The Capsid Proteins of a Large, Icosahedral dsDNA Virus

Chilo iridescent virus (CIV) is a large (∼ 1850 Å diameter) insect virus with an icosahedral, T = 147 capsid, a double-stranded DNA (dsDNA) genome, and an internal lipid membrane. The structure of CIV was determined to 13 Å resolution by means of cryoelectron microscopy (cryoEM) and three-dimensional image reconstruction. A homology model of P50, the CIV major capsid protein (MCP), was built based on its amino acid sequence and the structure of the homologous Paramecium bursaria chlorella virus 1 Vp54 MCP. This model was fitted into the cryoEM density for each of the 25 trimeric CIV capsomers per icosahedral asymmetric unit. A difference map, in which the fitted CIV MCP capsomers were subtracted from the CIV cryoEM reconstruction, showed that there are at least three different types of minor capsid proteins associated with the capsomers outside the lipid membrane. “Finger” proteins are situated at many, but not all, of the spaces between three adjacent capsomers within each trisymmetron, and “zip” proteins are situated between sets of three adjacent capsomers at the boundary between neighboring trisymmetrons and pentasymmetrons. Based on the results of segmentation and density correlations, there are at least eight finger proteins and three dimeric and two monomeric zip proteins in one asymmetric unit of the CIV capsid. These minor proteins appear to stabilize the virus by acting as intercapsomer cross-links. One transmembrane “anchor” protein per icosahedral asymmetric unit, which extends from beneath one of the capsomers in the pentasymmetron to the internal leaflet of the lipid membrane, may provide additional stabilization for the capsid. These results are consistent with the observations for other large, icosahedral dsDNA viruses that also utilize minor capsid proteins for stabilization and for determining their assembly.     

EADock: docking of small molecules into protein active sites with a multiobjective evolutionary optimization

In recent years, protein-ligand docking has become a powerful tool for drug development. Although several approaches suitable for high throughput screening are available, there is a need for methods able to identify binding modes with high accuracy. This accuracy is essential to reliably compute the binding free energy of the ligand. Such methods are needed when the binding mode of lead compounds is not determined experimentally but is needed for structure-based lead optimization. We present here a new docking software, called EADock, that aims at this goal. It uses an hybrid evolutionary algorithm with two fitness functions, in combination with a sophisticated management of the diversity. EADock is interfaced with the CHARMM package for energy calculations and coordinate handling. A validation was carried out on 37 crystallized protein-ligand complexes featuring 11 different proteins. The search space was defined as a sphere of 15 A around the center of mass of the ligand position in the crystal structure, and on the contrary to other benchmarks, our algorithm was fed with optimized ligand positions up to 10 A root mean square deviation (RMSD) from the crystal structure, excluding the latter. This validation illustrates the efficiency of our sampling strategy, as correct binding modes, defined by a RMSD to the crystal structure lower than 2 A, were identified and ranked first for 68% of the complexes. The success rate increases to 78% when considering the five best ranked clusters, and 92% when all clusters present in the last generation are taken into account. Most failures could be explained by the presence of crystal contacts in the experimental structure. Finally, the ability of EADock to accurately predict binding modes on a real application was illustrated by the successful docking of the RGD cyclic pentapeptide on the alphaVbeta3 integrin, starting far away from the binding pocket.     

Structural basis of poly(3-hydroxybutyrate) hydrolysis by PhaZ7 depolymerase from Paucimonas lemoignei

The crystal structure of poly(3-hydroxybutyrate) (PHB) depolymerase PhaZ7 purified from Paucimonas lemoignei was determined at 1.90 Å resolution. The structure consists of a single domain with an α/β hydrolase fold in its core. The active site is analogous to that of serine esterases/lipases and is characterized by the presence of a catalytic triad comprising Ser136, Asp242, and His306. Comparison with other structures in the Protein Data Bank showed a high level of similarity with the Bacillus subtilis lipase LipA (RMSD, 1.55 Å). Structural comparison with Penicillium funiculosum PHB depolymerase, the only PHB depolymerase whose structure is already known, revealed significant differences, resulting in an RMSD of 2.80-3.58 Å. The two enzymes appear to utilize different types of solvent-exposed residues for biopolymer binding, with aliphatic and hydroxyl residues used in P. funiculosum PHB depolymerase and aromatic residues in PhaZ7. Moreover, the active site of P. funiculosum PHB depolymerase is accessible to the substrate in contrast to the active site of PhaZ7, which is buried. Hence, considerable conformational changes are required in PhaZ7 for the creation of a channel leading to the active site. Taken together, the structural data suggest that PhaZ7 and P. funiculosum PHB depolymerase have adopted different strategies for effective substrate binding in response to their diverse substrate specificity and the lack of a substrate-binding domain.     

Comparing protein-ligand docking programs is difficult

There is currently great interest in comparing protein-ligand docking programs. A review of recent comparisons shows that it is difficult to draw conclusions of general applicability. Statistical hypothesis testing is required to ensure that differences in pose-prediction success rates and enrichment rates are significant. Numerical measures such as root-mean-square deviation need careful interpretation and may profitably be supplemented by interaction-based measures and visual inspection of dockings. Test sets must be of appropriate diversity and of good experimental reliability. The effects of crystal-packing interactions may be important. The method used for generating starting ligand geometries and positions may have an appreciable effect on docking results. For fair comparison, programs must be given search problems of equal complexity (e.g. binding-site regions of the same size) and approximately equal time in which to solve them. Comparisons based on rescoring require local optimization of the ligand in the space of the new objective function. Re-implementations of published scoring functions may give significantly different results from the originals. Ostensibly minor details in methodology may have a profound influence on headline success rates.     

Structure of the Small Outer Capsid Protein, Soc: A Clamp for Stabilizing Capsids of T4-like Phages

Many viruses need to stabilize their capsid structure against DNA pressure and for survival in hostile environments. The 9-kDa outer capsid protein (Soc) of bacteriophage T4, which stabilizes the virus, attaches to the capsid during the final stage of maturation. There are 870 Soc molecules that act as a “glue” between neighboring hexameric capsomers, forming a “cage” that stabilizes the T4 capsid against extremes of pH and temperature. Here we report a 1.9 Å resolution crystal structure of Soc from the bacteriophage RB69, a close relative of T4. The RB69 crystal structure and a homology model of T4 Soc were fitted into the cryoelectron microscopy reconstruction of the T4 capsid. This established the region of Soc that interacts with the major capsid protein and suggested a mechanism, verified by extensive mutational and biochemical studies, for stabilization of the capsid in which the Soc trimers act as clamps between neighboring capsomers. The results demonstrate the factors involved in stabilizing not only the capsids of T4-like bacteriophages but also many other virus capsids.     

Crystallographic Snapshot of Glycosylasparaginase Precursor Poised for Autoprocessing

Glycosylasparaginase belongs to a family of N-terminal nucleophile hydrolases that autoproteolytically generate their mature enzymes from single-chain protein precursors. Previously, based on a precursor structure paused at pre-autoproteolysis stage by a reversible inhibitor (glycine), we proposed a mechanism of intramolecular autoproteolysis. A key structural feature, a highly strained conformation at the scissile peptide bond, had been identified and was hypothesized to be critical for driving autoproteolysis through an N-O acyl shift. To examine this “twist-and-break” hypothesis, we report here a 1. 9-Å-resolution structure of an autoproteolysis-active precursor (a T152C mutant) that is free of inhibitor or ligand and is poised to undergo autoproteolysis. The current crystallographic study has provided direct evidence for the natural conformation of the glycosylasparaginase autocatalytic site without influence from any inhibitor or ligand. This finding has confirmed our previous proposal that conformational strain is an intrinsic feature of an active precursor.     

Crystal Structures of Limulus SAP-Like Pentraxin Reveal Two Molecular Aggregations

The serum-amyloid-P-component-like pentraxin from Limulus polyphemus, a recently discovered pentraxin species and important effector protein of the hemolymph immune system, displays two distinct doubly stacked cyclic molecular aggregations, heptameric and octameric. The refined three-dimensional structures determined by X-ray crystallography, both based on the same cDNA sequence, show that each aggregate is constructed from a similar dimer of protomers, which is repeated to make up the ring structure. The native octameric form has been refined at a resolution of 3 Å, the native heptameric form at 2.3 Å, and the phosphoethanolamine (PE)-bound octameric form at 2.7 Å. The existence of the hitherto undescribed heptameric form was confirmed by single-particle analysis using cryo-electron microscopy. In the native structures, the calcium-binding site is similar to that in human pentraxins, with two calcium ions bound in each subunit. Upon binding PE, however, each subunit binds a third calcium ion, with all three calcium ions contributing to the binding and orientation of the bound phosphate group within the ligand-binding pocket. While the phosphate is well-defined in the electron density, the ethanolamine group is poorly defined, suggesting structural and binding variabilities of this group. Although sequence homology with human serum amyloid P component is relatively low, structural homology is high, with very similar overall folds and a common affinity for PE. This is due, in part, to a “topological” equivalence of side-chain position. Identical side chains that are important in both function and fold, from different regions of the sequence in human and Limulus structures, occupy similar space within the overall subunit fold. Sequence and structure alignment, based on the refined three-dimensional structures presented here and the known horseshoe crab pentraxin sequences, suggest that adaptation and refinement of C-reactive-protein-mediated immune responses in these ancient creatures lacking antibody-based immunity are based on adaptation by gene duplication.     

Energy-Based Analysis and Prediction of the Orientation between Light- and Heavy-Chain Antibody Variable Domains

Diversity in antibody structure is crucial to the ability of the adaptive immune system to recognize the tremendously diverse set of potential antigens. The diversity in structure is most apparent in the six hypervariable loops of the complementarity-determining regions. However, given that these loops occur at the interface of the heavy- and light-chain variable domains and form the antigen-binding site, the relative orientation of the heavy- and light-chain variable domains can create another source of structural diversity leading to changes in antigen binding. Here, we first reexamine the diversity of V_L:V_H orientations in existing antibody crystal structures using 153 nonredundant sequences, demonstrating that the variation in V_L:V_H orientation is greater than that expected from effects of crystal packing, antigen binding, or the presence of antibody constant regions and increases, on average, as sequence similarity decreases for residues in the interface between the domains. We developed a tool for predicting the relative orientations of the heavy- and light-chain variable domains using side-chain rotamer sampling in the interface and molecular-mechanics-based energy calculations. When using variable domain backbones from the crystal structures, the predicted orientation is very close (< 1 Å RMSD) to the crystallographically observed orientation in most cases, confirming that the V_L:V_H orientation is determined by the antibody sequence and suggesting an approach to predicting the relative orientation of the variable domains when building homology models of antibodies. When applied to antibody homology models generated from templates with 55-75% sequence identity, we predict the V_L:V_H orientation of 20 antibodies with an average/median RMSD of 2.1/1.6 Å to the crystal structures.     

Antiparallel Conformation of Knob and Hole Aglycosylated Half-Antibody Homodimers Is Mediated by a CH2–CH3 Hydrophobic Interaction

Bispecific antibody and antibody-like molecules are of wide interest as potential therapeutics that can recognize two distinct targets. Among the variety of ways such molecules have been engineered is by creating “knob” and “hole” heterodimerization sites in the CH3 domains of two antibody heavy chains. The molecules produced in this manner maintain their biological activities while differing very little from the native human IgG sequence. To better understand the knob-into-hole interface, the molecular mechanism of heterodimerization, and to engineer Fc domains that could improve the assembly and purity of heterodimeric reaction products, we sought crystal structures of aglycosylated heterodimeric and homodimeric “knob” and “hole” Fc fragments derived from bacterial expression. The structure of the knob-into-hole Fc was determined at 2.64 Å. Except for the sites of mutation, the structure is very similar to that of the native human IgG1 Fc, consistent with a heterodimer interaction kinetic K_D of < 1 nM. Homodimers of the “knob” and “hole” mutants were also obtained, and their X-ray structures were determined at resolutions 2.5 Å and 2.1 Å, respectively. Both kinds of homodimers adopt a head-to-tail quaternary structure and thus do not contain direct knob/knob or hole/hole CH3 interactions. The head-to-tail arrangement was disfavored by adding site-directed mutations at F241 and F243 in the CH2 domains, leading to increases in both rate and efficiency of bispecific (heterodimer) assembly.     

Cytoplasmic oxysterol-binding proteins: sterol sensors or transporters?

Families of oxysterol-binding protein (ORP) homologues are present in eukaryotes from yeast to man. Their hallmark feature is a characteristic ligand binding domain that, for several family members, has been shown to accommodate different oxysterols and/or cholesterol. ORPs of the “long” subtype contain targeting determinants for the endoplasmic reticulum and to other organelle membranes, the most prominent of which are phosphoinositide-binding pleckstrin homology domains, while “short” ORPs comprise a ligand binding domain with little additional sequences. There is increasing evidence that both long and short ORPs can be enriched at membrane contact sites, junctions of the endoplasmic reticulum with other organelles, where they are suggested to execute regulatory or sterol transfer functions. In this review we discuss the current evidence for putative roles of ORPs as sterol sensors or transporters.     

ConsDock: A new program for the consensus analysis of protein-ligand interactions

Protein-based virtual screening of chemical libraries is a powerful technique for identifying new molecules that may interact with a macromolecular target of interest. Because of docking and scoring limitations, it is more difficult to apply as a lead optimization method because it requires that the docking/scoring tool is able to propose as few solutions as possible and all of them with a very good accuracy for both the protein-bound orientation and the conformation of the ligand. In the present study, we present a consensus docking approach (ConsDock) that takes advantage of three widely used docking tools (Dock, FlexX, and Gold). The consensus analysis of all possible poses generated by several docking tools is performed sequentially in four steps: (i) hierarchical clustering of all poses generated by a docking tool into families represented by a leader; (ii) definition of all consensus pairs from leaders generated by different docking programs; (iii) clustering of consensus pairs into classes, represented by a mean structure; and (iv) ranking the different means starting from the most populated class of consensus pairs. When applied to a test set of 100 protein-ligand complexes from the Protein Data Bank, ConsDock significantly outperforms single docking with respect to the docking accuracy of the top-ranked pose. In 60% of the cases investigated here, ConsDock was able to rank as top solution a pose within 2 A RMSD of the X-ray structure. It can be applied as a postprocessing filter to either single- or multiple-docking programs to prioritize three-dimensional guided lead optimization from the most likely docking solution.