Prediction of binding constants of protein ligands: A fast method for the prioritization of hits obtained from de novo design or 3D database search programs

A dataset of 82 protein–ligand complexes of known 3D structure and binding constant K_i was analysed to elucidate the important factors that determine the strength of protein–ligand interactions. The following parameters were investigated: the number and geometry of hydrogen bonds and ionic interactions between the protein and the ligand, the size of the lipophilic contact surface, the flexibility of the ligand, the electrostatic potential in the binding site, water molecules in the binding site, cavities along the protein–ligand interface and specific interactions between aromatic rings. Based on these parameters, a new empirical scoring function is presented that estimates the free energy of binding for a protein–ligand complex of known 3D structure. The function distinguishes between buried and solvent accessible hydrogen bonds. It tolerates deviations in the hydrogen bond geometry of up to 0.25 Å in the length and up to 30 °Cs in the hydrogen bond angle without penalizing the score. The new energy function reproduces the binding constants (ranging from 3.7 × 10^-2 M to 1 × 10^-14 M, corresponding to binding energies between -8 and -80 kJ/mol) of the dataset with a standard deviation of 7.3 kJ/mol corresponding to 1.3 orders of magnitude in binding affinity. The function can be evaluated very fast and is therefore also suitable for the application in a 3D database search or de novo ligand design program such as LUDI. The physical significance of the individual contributions is discussed.     

Properties of Sulfur-Sulfur Bonds

SS bonds are extraordinarily flexible and have properties that are observed only on isolated occasions for other homonuclear bonds: the bond lengths very between 1.8 and 3.0Å, the bond angles between 90 and 180° and the dihedral angles between 0 and 180°; the bond energies amount to up to 430 kJ/mol. The SS stretching frequencies can appear over the range 177–820 cm^?1 and force constants of 1.4 to 6.3 mdyne/Å have been calculated. This variability is illustrated with examples containing isolated and cumulated SS bonds.     

Torsion angle preference and energetics of small-molecule ligands bound to proteins

Small organic molecules can assume conformations in the protein-bound state that are significantly different from those in solution. We have analyzed the conformations of 21 common torsion motifs of small molecules extracted from crystal structures of protein-ligand complexes and compared them with their torsion potentials calculated by an ab initio DFT method. We find a good correlation between the potential energy of the torsion motifs and their conformational distribution in the protein-bound state: The most probable conformations of the torsion motifs agree well with the calculated global energy minima, and the lowest torsion-energy state becomes increasingly dominant as the torsion barrier height increases. The torsion motifs can be divided into 3 groups based on torsion barrier heights: high (>4 kcal/mol), medium (2-4 kcal/mol), and low (<2 kcal/mol). The calculated torsion energy profiles are predictive for the most preferred bound conformation for the high and medium barrier groups, the latter group common in druglike molecules. In the high-barrier group of druglike ligands, >95% of conformational torsions occur in the energy region <4 kcal/mol. The conformations of the torsion motifs in the protein-bound state can be modeled by a Boltzmann distribution with a temperature factor much higher than room temperature. This high-temperature factor, derived by fitting the theoretical model to the experimentally observed conformation occurrence of torsions, can be interpreted as the perturbation that proteins inflict on the conformation of the bound ligand. Using this model, it is calculated that the average strain energy of a torsion motif in ligands bound to proteins is approximately 0.6 kcal/mol, a result which can be related to the lower binding efficiency of larger ligands with more rotatable bonds. The above results indicate that torsion potentials play an important role in dictating ligand conformations in both the free and the bound states.     

Stabilities of Hydrogen-Bonded Supramolecular Complexes with Various Numbers of Single Bonds: Attempts To Quantify a Dogma in Host–Guest Chemistry

The disadvantage of flexible bonds in supramolecular host–guest complexes is much smaller than usually assumed. In the association of dicarboxlic acids with diamides (shown on the right), freely rotatable single bonds have only a minor disadvantageous influence on the free energy of complexation ΔG_cplx. From the linear correlation between ΔG_cplx and the number of single bonds n, a decrease in the free energy of complexation of only 1.3 kJ mol⁻¹ per single bond was determined.     

Cocrystal Assembled by 1,2-Diiodotetrafluorobenzene and Acridine via C-I…N Halogen Bond and Tr-hole…F Bonds

The concepts on o-hole and ~-hole bonds are suggested. A cocrystal with repeated 8-F-atom unit as basic struc- tural motif is assembled based on bifurcated C-I…N…I-C halogen/σ-hole bond and antiparallel double π-hole… F bonds by 1,2-diiodotetrafluorobenzene and acridine and characterized well by XRD, powder XRD and solid 19F NMR, etc. Also the calculated interaction energies are -26.8 and -31.5 kJ/mol for bifurcated C-I…N sp……2 halogen bonds, and -14.3 kJ/mol for a pair of n-hole…F bonds. In this system C-I…N halogen bond has stronger competitive ability to C-I…π halogen bond due to stronger basicity of N than π-system in acridine. The combination of the halogen/σ-hole and π-hole bonds or together with other weak interactions could play a key role in assembling function materials, molecular recognition and design of drugs and so on.     

Quantum-Chemical Study of Adducts of Silicon Halides with Nitrogen-containing Donors: IV. Adducts with Pyridine

The structural and thermodynamic characteristics of SiX₄·Py and SiX₄·2Py adducts (X = H, F, Cl, Br) were calculated by ab initio and DFT methods (RHF and B3LYP). The resulting data were used to estimate for the first time the enthalpies of sublimation of trans-SiX₄·2Py complexes. The distortion energies of the donor and acceptor fragments and the energies of the Si-N bonds in the 1:1 and 1:2 halide complexes were calculated. The high distortion energy makes thermodynamically unfavorable equatorial monopyridine adducts with Si-N bond energies of 150-200 kJ/mol. In trans 1:2 complexes, pyridine acts as a weaker donor than ammonia with respect to silicon tetrahalides.     

Computational evidence for methyl-donated hydrogen bonds and hydrogen-bond networking in 1,2-ethanediol-dimethyl sulfoxide

The 1:1 complex of 1,2-ethanediol with dimethyl sulfoxide was studied using density functional theory. A network of three hydrogen bonds holds the complex together, including two in which each methyl group donates to the same hydroxyl oxygen. Four lines of evidence support the existence of methyl-donated hydrogen bonds. The interaction energy is 36 +/- 5 kJ/mol using Becke's three parameter hybrid theory with the 1991 nonlocal correlation functional of Perdew and Wang, and a moderately large basis set (B3PW91/6-311++G**//B3PW91/6-31+G**). To determine the energy of each hydrogen bond, a relaxed potential energy scan was performed in a smaller basis set to break the weaker hydrogen bonds by forced systematic rotation of the methyl groups. Two cross-checking analyses show cooperative effects that cause individual hydrogen bond energies in the network to be nonadditive. When one methyl hydrogen bond is broken, the remaining interactions stabilize the complex by storing an additional 2-3 kJ/mol. With all hydrogen bonds intact, the O[bond]H...O[bond]S hydrogen bond contributes 26 +/- 2 kJ/mol stability, and each weak methyl bond stores 5 +/- 2 kJ/mol.     

Kinetics of the Dissolution of Silver in Thiocarbamide Solutions

The dependence of the rate of solution of silver on the pH of the solution, the ratio of the iron(III) and thiocarbamide concentrations, and the temperature has been determined. The rate constants for the solution of silver (k _i = 2.3·10^–4 to 9.6·10^–4s^–1) at temperatures from 283-298 K have been calculated and from the temperature dependence of the rate constant the activation energies have been calculated: 68.84 kJ/mol for kinetic control of the rate of solution and 26.06 kJ/mol in the adsorption inhibition region.     

Formation of MSeO4(aq) complexes (M2+ = Mg2+, Co2+, Ni2+, Cu2+, Cd2+) studied as a function of temperature by affinity capillary electrophoresis

Complexation of divalent cations (Mg²⁺, Co²⁺, Ni²⁺, Cu²⁺, Cd²⁺) by selenate ligand was studied by ACE (UV indirect detection) in 0.1 mol/L NaNO₃ ionic strength solutions at various temperatures (15, 25, 35, 45 and 55°C). For each solution, a unique peak was observed as a result of a fast equilibrium between the free ion and the complex (labile systems). The migration time corresponding to this peak changed as a function of the solution composition, namely the free and complexed metal concentrations, according to the complexation reactions. The results confirmed the formation of a unique 1:1 complex for each cation. The thermodynamic parameters were fitted to the experimental data at 0.1 mol/L ionic strength: (25°C) = ?(6.5 ± 0.3), ?(7.5 ± 0.3), ?(7.7 ± 0.3), ?(7.7 ± 0.3), and –(8.1 ± 0.3) kJ/mol and = 2.5 ± 0.2, 4.7 ± 0.4, 4.5 ± 0.6, 8.4 ± 1.1, and 7.2 ± 0.6 kJ/mol for M²⁺ = Mg²⁺, Co²⁺, Ni²⁺, Cu²⁺, and Cd²⁺, respectively. Complexes with alkaline earth and transition metal cations could be distinguished by their relative stabilities. The effect of the ionic medium was treated using the specific ion interaction theory and the thermodynamic parameters at infinite dilution were compared to previously published data on metal–selenate, metal–sulfate, and metal–chromate complexes.     

Thermogravimetric and differential scanning calorimetric investigations of manganese–glycine interactions

Thermogravimetric (t.g.) and differential scanning calorimetric (d.s.c.) data have been used to study metal–amino acid interactions in adducts of general formula MnCl₂ · ngly (gly = glycine, n = 0.7, 2.0, 4.0 and 5.0). All the prepared adducts exhibit only a one step mass loss associated with the release of glycine molecules, except for the 0.7gly adduct, which exhibits two glycine mass loss steps. From d.s.c. data, the enthalpy values associated with the glycine mass loss can be calculated: MnCl₂ · 0.7gly = 409 and 399 kJ mol^–1, MnCl₂ · 2.0gly = 216 kJ mol^–1, MnCl₂ · 4.0gly = 326 kJ mol^–1 and MnCl₂ · 5.0gly = 423 kJ mol^–1, respectively. The enthalpy associated with the ligand loss, plotted as function of the number of ligands for the n = 2.0, 4.0 and 5.0 adducts, gave a linear correlation, fitting the equation: H (ligand loss)/kJ mol^–1 = 67 × (number of ligands, n) + 76. A similar result was achieved when the enthalpy associated with the ligand loss was plotted as a function of the _a(COO^–) bands associated with the coordination through the carboxylate group, 1571, 1575 and 1577 cm^–1, respectively, for the n = 2.0, 4.0 and 5.0 adducts, giving the equation H (ligand loss) /kJ mol^–1 = 33.5 × _a(COO^–) /cm^–1 – 52418.5. This simple equation provides evidence for the enthalpy associated with the ligand loss being very closely related to the electronic density associated with the metal–amino acid bonds.     

Molecular structure and binding energies of monosubstituted hexacarbonyls of chromium,molybdenum, and tungsten: Relativistic density functional study

Relativistic density functional calculations have been carried out for the group VI transition metal carbonyls M(CO)₅L (M=Cr, Mo, W; L=OH₂, NH₃, PH₃, PMe₃, N₂, CO, OC (isocarbonyl), CS, CH₂, CF₂, CCl₂, NO⁺). The optimized molecular structures and M(SINGLE BOND)L bond dissociation energies, as well as the metal–carbonyl bond energy of the trans CO group, have been calculated. Besides the marked dependence of the trans M(SINGLE BOND)CO bond length on the type of ligand L, such an effect on the that bond energy is also observed. For the chromium compounds, the trans Cr(SINGLE BOND)CO bond length varies from 184 to 199 pm and its bond energy from 242 to 150 kJ/mol. For the molybdenum compounds, the range is 197 to 216 pm and 253 to 128 kJ/mol and, for tungsten, 198 to 214 pm and 293 to 159 kJ/mol. The observed trends can be explained with the π acceptor strength of the L ligand. © 1997 John Wiley & Sons, Inc. J Comput Chem 18 : 1985–1992, 1997     

A Mixed QM/MM Scoring Function to Predict Protein-Ligand Binding Affinity

Computational methods for predicting protein-ligand binding free energy continue to be popular as a potential cost-cutting method in the drug discovery process. However, accurate predictions are often difficult to make as estimates must be made for certain electronic and entropic terms in conventional force field based scoring functions. Mixed quantum mechanics/molecular mechanics (QM/MM) methods allow electronic effects for a small region of the protein to be calculated, treating the remaining atoms as a fixed charge background for the active site. Such a semi-empirical QM/MM scoring function has been implemented in AMBER using DivCon and tested on a set of 23 metalloprotein-ligand complexes, where QM/MM methods provide a particular advantage in the modeling of the metal ion. The binding affinity of this set of proteins can be calculated with an R(2) of 0.64 and a standard deviation of 1.88 kcal/mol without fitting and 0.71 and a standard deviation of 1.69 kcal/mol with fitted weighting of the individual scoring terms. In this study we explore using various methods to calculate terms in the binding free energy equation, including entropy estimates and minimization standards. From these studies we found that using the rotational bond estimate to ligand entropy results in a reasonable R(2) of 0.63 without fitting. We also found that using the ESCF energy of the proteins without minimization resulted in an R(2) of 0.57, when using the rotatable bond entropy estimate.     

Conformers and intramolecular hydrogen bonding of the oxalic acid monomer and its anions

Density functional theory, B3LYP/6‐31G** and B3LYP/6‐311+G(2d,p), and ab initio MP2/6‐31G** calculations have been carried out to investigate the conformers, transition states, and energy barriers of the conformational processes of oxalic acid and its anions. QCISD/6‐31G** geometrical optimization is also performed in the stable forms. Its calculated energy differences between the two most stable conformers are very near to the related observed value at 7.0 kJ/mol. It is found that the structures and relative energies of oxalic acid conformers predicted by these methods show similar results, and that the conformer L1 (C_2h) with the double‐interfunctional‐groups hydrogen bonds is the most stable conformer. The magnitude of hydrogen bond energies depends on the energy differences of various optimized structures. The hydrogen bond energies will be about 32 kJ/mol for interfunctional groups, 17 kJ/mol for weak interfunctional groups, 24 kJ/mol for intra‐COOH in (COOH)₂, and 60 kJ/mol for interfunctional groups in (COOH)COO⁻¹ ion if calculated using the B3LYP/6‐311+G(2d,p) method. © 2000 John Wiley & Sons, Inc. Int J Quant Chem 76: 541–551, 2000     

Hydration in drug design. 1. Multiple hydrogen-bonding features of water molecules in mediating protein-ligand interactions

Summary Water is known to play an important rôle in the recognition and stabilization of the interaction between a ligand and its site. This has important implications for drug design. Analyses of 19 high-resolution crystal structures of protein-ligand complexes reveal the multiple hydrogen-bonding feature of water molecules mediating protein-ligand interactions. Most of the water molecules (nearly 80%) involved in bridging the protein and the ligand can make three or more hydrogen bonds when distance and bond angles are used as criteria to define hydrogen-bonding interactions. Isotropic B-factors have been used to take into account the mobility of water molecules. The water molecules at binding sites bridge the protein and ligand, and interact with other water molecules to form a complex network of interconnecting hydrogen bonds. Some water molecules at the site do not directly bridge between the protein and the ligand, but may contribute indirectly to the stability of the complex by holding bridging water molecules in the right position through a network of hydrogen bonds. These water networks are probably crucial for the stability of the protein-ligand complex and are important for any site-directed drug design strategies.     

Binding affinities for models of biologically available potential Cu(II) ligands relevant to Alzheimer's disease: an ab initio study

A systematic study of the binding affinities of the model biological ligands X: = (CH3)2S, CH3S-, CH3NH2, 4-CH3-imidazole (MeImid), C6H5O-, and CH3CO2- to (NH3)i(H2O)3-iCu(II)-H2O (i = 3, 2, 1, 0) complexes has been carried out using quantum chemical calculations. Geometries have been obtained at the B3LYP/ 6-31G(d) level of theory, and binding energies, Delta, relative to H2O as a ligand, have been calculated at the B3LYP/6-311+G(2df,2p)//B3LYP/6-31G(d) level. Solvation effects have been included using the COSMO model, and the relative binding free energies in aqueous solution (Delta) have been determined at pH 7 for processes that are pH dependent. CH3S- (Delta = -16.0 to -53.5 kJ mol(-1)) and MeImid (Delta = -18.5 to -35.2 kJ mol(-1)) give the largest binding affinities for Cu(II). PhO- and (CH3)2S are poor ligands for Cu(II), Delta = 20.6 to -9.7 and 19.8 to -3.7 kJ mol(-1), respectively. The binding affinities for CH3NH2 range from -0.8 to -15.0 kJ mol(-1). CH3CO2- has Cu(II) binding affinities in the ranges Delta = -13.5 to -32.4 kJ mol(-1) if an adjacent OH bond is available for hydrogen bonding and Delta = 10.1 to -4.6 kJ mol(-1) if this interaction is not present. In the context of copper coordination by the Abeta peptide of Alzheimer's disease, the binding affinities suggest preferential binding of Cu(II) to the three histidine residues plus a lysine or the N-terminus. For a 3N1O Cu(II) ligand arrangement, it is more probable that the oxygen ligand comes from an aspartate/glutamate residue side chain than from the tyrosine at position 10. Methionine appears unlikely to be a Cu(II) ligand in Abeta.     

Calculations of protein-ligand binding entropy of relative and overall molecular motions

In the context of virtual database screening, calculations of protein-ligand binding entropy of relative and overall molecular motions are challenging, owing to the inherent structural complexity of the ligand binding well in the energy landscape of protein-ligand interactions and computing time limitations. We describe a fast statistical thermodynamic method for estimation the binding entropy to address the challenges. The method is based on the integration of the configurational integral over clusters obtained from multiple docked positions. We apply the method in conjunction with 11 popular scoring functions (AutoDock, ChemScore, DrugScore, D-Score, F-Score, G-Score, LigScore, LUDI, PLP, PMF, X-Score) to evaluate the binding entropy of 100 protein-ligand complexes. The averaged values of binding entropy contribution vary from 6.2 to 9.1 kcal/mol, showing good agreement with literature. We calculate positional sizes and the angular volume of the native ligand wells. The averaged geometric mean of positional sizes in principal directions varies from 0.8 to 1.4 A. The calculated range of angular volumes is 3.3-11.8 rad(2). Then we demonstrate that the averaged six-dimensional volume of the native well is larger than the volume of the most populated non-native well in energy landscapes described by all of 11 scoring functions.     

Theoretical study of the conformational states of triosmium clusters with a chiral µ-1-NH pinane ligand

This paper reports a theoretical investigation of the conformations of triosmium clusters with a chiral pinane ligand (μ-H)Os₃(CO)₁₀(μ-NHC₁₀H₁₇). The potential curves of internal rotation of the organic ligand relative to the N-C bond are plotted for the cluster complexes. The structures of possible conformers are considered, and reasons for their stability are revealed. The barrier of rotation around the N-C bond of the terpenoid is 186.6 kJ/mol for crystals and ∼140 kJ/mol for solutions. Due to this, the free rotation of the ligand is hindered in both cases. The effects of the intra-and intermolecular interactions on the conformational state of the cluster complex are analyzed. Original Russian Text Copyright ? 2007 by V. A. Potemkin, V. A. Maksakov, and V. S. Korenev __________ Translated from Zhurnal Strukturnoi Khimii, Vol. 48, No. 2, pp. 230–235, March–April, 2007.     

The nature of interactions between N-butylpyridinium tetrafluoroborate and thiophenic compounds: A theoretical investigation

In an effort to understand the nature of the interactions between pyridinium-based ionic liquids and thiophenic compounds, the electronic and topological properties of the interactions between N-butylpyridinium tetrafluoroborate ([BPY]⁺[BF₄]⁻) and thiophene (TS), benzothiophene (BT), dibenzothiophene (DBT) have been investigated by density functional theory. The most stable structure of the [BPY]⁺[BF₄]⁻ ion-pair indicated that hydrogen bonding interactions between fluorine atoms on [BF₄]⁻ anions and C2–H2 on the pyridinium ring play an important role in the formation of the ion-pair. The NBO and AIM analyses indicate the occurrence of π–π stacking interactions. The electron density at bond critical points and Wiberg bond indices are correlated with the interacting distances of H···F interactions, so electron density and Wiberg bond index can demonstrate the interacting strength of H···F hydrogen bonds. The interaction energies suggest that DBT adsorbs prior to the other compounds on N-butylpyridinium tetrafluoroborate ionic liquid.     

Accelerating molecular docking calculations using graphics processing units

The generation of molecular conformations and the evaluation of interaction potentials are common tasks in molecular modeling applications, particularly in protein-ligand or protein-protein docking programs. In this work, we present a GPU-accelerated approach capable of speeding up these tasks considerably. For the evaluation of interaction potentials in the context of rigid protein-protein docking, the GPU-accelerated approach reached speedup factors of up to over 50 compared to an optimized CPU-based implementation. Treating the ligand and donor groups in the protein binding site as flexible, speedup factors of up to 16 can be observed in the evaluation of protein-ligand interaction potentials. Additionally, we introduce a parallel version of our protein-ligand docking algorithm PLANTS that can take advantage of this GPU-accelerated scoring function evaluation. We compared the GPU-accelerated parallel version to the same algorithm running on the CPU and also to the highly optimized sequential CPU-based version. In terms of dependence of the ligand size and the number of rotatable bonds, speedup factors of up to 10 and 7, respectively, can be observed. Finally, a fitness landscape analysis in the context of rigid protein-protein docking was performed. Using a systematic grid-based search methodology, the GPU-accelerated version outperformed the CPU-based version with speedup factors of up to 60.     

Empirical scoring functions: I. The development of a fast empirical scoring function to estimate the binding affinity of ligands in receptor complexes

This paper describes the development of a simple empirical scoringfunction designed to estimate the free energy of binding for aprotein–ligand complex when the 3D structure of the complex is knownor can be approximated. The function uses simple contact terms to estimatelipophilic and metal–ligand binding contributions, a simple explicitform for hydrogen bonds and a term which penalises flexibility. Thecoefficients of each term are obtained using a regression based on 82ligand–receptor complexes for which the binding affinity is known. Thefunction reproduces the binding affinity of the complexes with across-validated error of 8.68 kJ/mol. Tests on internal consistency indicatethat the coefficients obtained are stable to changes in the composition ofthe training set. The function is also tested on two test sets containing afurther 20 and 10 complexes, respectively. The deficiencies of this type offunction are discussed and it is compared to approaches by other workers.