Empirical free energy calculation: comparison to calorimetric data.

An effective free energy potential, developed originally for binding free energy calculation, is compared to calorimetric data on protein unfolding, described by a linear combination of changes in polar and nonpolar surface areas. The potential consists of a molecular mechanics energy term calculated for a reference medium (vapor or nonpolar liquid), and empirical terms representing solvation and entropic effects. It is shown that, under suitable conditions, the free energy function agrees well with the calorimetric expression. An additional result of the comparison is an independent estimate of the side-chain entropy loss, which is shown to agree with a structure-based entropy scale. These findings confirm that simple functions can be used to estimate the free energy change in complex systems, and that a binding free energy evaluation model can describe the thermodynamics of protein unfolding correctly. Furthermore, it is shown that folding and binding leave the sum of solute-solute and solute-solvent van der Waals interactions nearly invariant and, due to this invariance, it may be advantageous to use a nonpolar liquid rather than vacuum as the reference medium.     

Free energy distributions in proteins.

Protein molecules in solution have a broad distribution of enthalpy states. A good approximation to the distribution function for enthalpy states can be calculated, using the maximum-entropy method, from the moments of the distribution that, in turn, are obtained from the experimental temperature dependence of the heat capacity. In the present paper, we show that the enthalpy probability distribution can then be formulated in terms of a free energy function that gives the free energy of the protein corresponding to a particular value of the enthalpy. By the location of the minima in this function, the free energy distribution graphically indicates the most probable values of the enthalpy for the protein. We find that the behavior of the free energy functions for proteins falls somewhere between two different cases: a two-state like function with two minima, the relative levels of the two states changing with temperature; and, a single-minimum function where the position of the minimum shifts to higher enthalpy values as the temperature is increased. We show that the temperature dependence of the free energy function can be expressed in terms of a central free energy distribution for a given, fixed temperature (which is most conveniently chosen as the temperature of the maximum in the heat capacity). The nature of this central free energy function for a given protein thus yields all of the thermodynamic behavior of that protein over the temperature range of the denaturation process.     

Insights into protein-protein binding by binding free energy calculation and free energy decomposition for the Ras-Raf and Ras-RalGDS complexes

Absolute binding free energy calculations and free energy decompositions are presented for the protein-protein complexes H-Ras/C-Raf1 and H-Ras/RalGDS. Ras is a central switch in the regulation of cell proliferation and differentiation. In our study, we investigate the capability of the molecular mechanics (MM)-generalized Born surface area (GBSA) approach to estimate absolute binding free energies for the protein-protein complexes. Averaging gas-phase energies, solvation free energies, and entropic contributions over snapshots extracted from trajectories of the unbound proteins and the complexes, calculated binding free energies (Ras-Raf: -15.0(+/-6.3)kcal mol(-1); Ras-RalGDS: -19.5(+/-5.9)kcal mol(-1)) are in fair agreement with experimentally determined values (-9.6 kcal mol(-1); -8.4 kcal mol(-1)), if appropriate ionic strength is taken into account. Structural determinants of the binding affinity of Ras-Raf and Ras-RalGDS are identified by means of free energy decomposition. For the first time, computationally inexpensive generalized Born (GB) calculations are applied in this context to partition solvation free energies along with gas-phase energies between residues of both binding partners. For selected residues, in addition, entropic contributions are estimated by classical statistical mechanics. Comparison of the decomposition results with experimentally determined binding free energy differences for alanine mutants of interface residues yielded correlations with r(2)=0.55 and 0.46 for Ras-Raf and Ras-RalGDS, respectively. Extension of the decomposition reveals residues as far apart as 25A from the binding epitope that can contribute significantly to binding free energy. These "hotspots" are found to show large atomic fluctuations in the unbound proteins, indicating that they reside in structurally less stable regions. Furthermore, hotspot residues experience a significantly larger-than-average decrease in local fluctuations upon complex formation. Finally, by calculating a pair-wise decomposition of interactions, interaction pathways originating in the binding epitope of Raf are found that protrude through the protein structure towards the loop L1. This explains the finding of a conformational change in this region upon complex formation with Ras, and it may trigger a larger structural change in Raf, which is considered to be necessary for activation of the effector by Ras.     

Heat capacity-independent determination of differential free energy of stability between structurally homologous proteins

Under the assumption of equivalent heat capacity values, the differential free energy of stability for a pair of proteins midway between their thermal unfolding transition temperatures is shown to be independent of DeltaC(p) up to its cubic term in DeltaT(m). For model calculations reflecting the nearly 30 degrees C difference in T(m) for the adenylate kinases from the arctic bacterium Bacillus globisporus and the thermophilic bacterium Geobacillus stearothermophilus, the resultant error in estimating DeltaDeltaG by the formula 0.5 [DeltaS(T(m1))(1)+DeltaS(T(m2)) (2)] DeltaT(m) is less than 1%. Combined with the analogous thermal unfolding data for the adenylate kinase from Escherichia coli, these three homologous proteins exhibit T(m) and DeltaS(T(m)) values consistent with differential entropy and enthalpy contributions of equal magnitude. When entropy-enthalpy compensation holds for the differential free energy of stability, the incremental changes in T(m) values are shown to be proportionate to the changes in free energy.     

Binding free energy calculation with QM/MM hybrid methods for Abl-Kinase inhibitor

We report a Quantum mechanics/Molecular Mechanics–Poisson-Boltzmann/ Surface Area (QM/MM-PB/SA) method to calculate the binding free energy of c-Abl human tyrosine kinase by combining the QM and MM principles where the ligand is treated quantum mechanically and the rest of the receptor by classical molecular mechanics. To study the role of entropy and the flexibility of the protein ligand complex in a solvated environment, molecular dynamics calculations are performed using a hybrid QM/MM approach. This work shows that the results of the QM/MM approach are strongly correlated with the binding affinity. The QM/MM interaction energy in our reported study confirms the importance of electronic and polarization contributions, which are often neglected in classical MM-PB/SA calculations. Moreover, a comparison of semi-empirical methods like DFTB-SCC, PM3, MNDO, MNDO-PDDG, and PDDG-PM3 is also performed. The results of the study show that the implementation of a DFTB-SCC semi-empirical Hamiltonian that is derived from DFT gives better results than other methods. We have performed such studies using the AMBER molecular dynamic package for the first time. The calculated binding free energy is also in agreement with the experimentally determined binding affinity for c-Abl tyrosine kinase complex with Imatinib.     

State-dependent ideal free distributions

Summary The standard ideal free distribution (IFD) states how animals should distribute themselves at a stable competitive equilibrium. The equilibrium is stable because no animal can increase its fitness by changing its location. In applying the IFD to choice between patches of food, fitness has been identified with the net rate of energetic gain. In this paper we assess fitness in terms of survival during a non-reproductive period, where the animal may die as a result of starvation or predation. We find the IFD when there is a large population that can distribute itself between two patches of food. The IFD in this case is state-dependent, so that an animal's choice of patch depends on its energy reserves. Animals switch between patches as their reserves change and so the resulting IFD is a dynamic equilibrium. We look at two cases. In one there is no predation and the patches differ in their variability. In the other, patches differ in their predation risk. In contrast to previous IFDs, it is not necessarily true that anything is equalized over the two patches.     

Exploring the interaction between D-xylose isomerase and D-xylose by free energy calculation

The numerical quadrature thermodynamic integration method is used to investigate enzyme-substrate interaction of D-xylose isomerase. A screening function for the coulombic interaction is introduced into the simulation to correct the effect of finite cutoff radius for the non-bonded interaction. The binding free energy difference for D-xylose with D-xylose isomerase and its N184D mutant has been calculated, and the result 3.9 ± 1.2 kJ/mol agrees well with experimental data of 4.38 kJ/mol. In addition, the structure and dynamics of enzyme-substrate complex were simulated for mutant and wild-type enzyme, respectively. Analysis of the structures and intramolecular interactions of the complexes were found to be valuable for understanding the reaction mechanism of the enzyme D-xylose isomerase. © Wiley-Liss, Inc.     

Improved free energy parameters for RNA pseudoknotted secondary structure prediction

Accurate prediction of RNA pseudoknotted secondary structures from the base sequence is a challenging computational problem. Since prediction algorithms rely on thermodynamic energy models to identify low-energy structures, prediction accuracy relies in large part on the quality of free energy change parameters. In this work, we use our earlier constraint generation and Boltzmann likelihood parameter estimation methods to obtain new energy parameters for two energy models for secondary structures with pseudoknots, namely, the Dirks–Pierce (DP) and the Cao–Chen (CC) models. To train our parameters, and also to test their accuracy, we create a large data set of both pseudoknotted and pseudoknot-free secondary structures. In addition to structural data our training data set also includes thermodynamic data, for which experimentally determined free energy changes are available for sequences and their reference structures. When incorporated into the HotKnots prediction algorithm, our new parameters result in significantly improved secondary structure prediction on our test data set. Specifically, the prediction accuracy when using our new parameters improves from 68% to 79% for the DP model, and from 70% to 77% for the CC model.     

Induced fit docking,free energy calculation and molecular dynamics studies on Mycobacterium tuberculosis alanine racemase inhibitor

Alar, a Pyridoxal 5′-phosphate (PLP)-dependent bacterial enzyme is responsible for the racemisation of L-alanine into D-alanine which is essential for the peptidoglycan biosynthesis in both Gram-positive and Gram-negative bacteria. In the present study, we performed induced fit docking, binding free energy calculation and molecular dynamics simulation to elucidate the Alar inhibition potential of 1,2,4-thiadiazolidine-3,5-dione-based inhibitor 1. The inhibitor binds to the hydrophobic groove of Alar and the binding was found to be stable throughout 20-ns MD simulation. Induced fit docking result showed that Lys42, Tyr46, Tyr175 and Tyr364 residues are primarily responsible for the stabilisation of inhibitor–protein complex. Further, high negative van der Waals binding free energy value of –38.88 kcal/mol, indicated it as the main driving force for the inhibitor binding. Based on the information obtained from this study, we designed few molecules as potent Alar inhibitor. In order to gain structural insight and to validate the stability of complex, we performed 20-ns MD simulation of the designed molecule D1. Results obtained from this study can be used for the design of M. tuberculosis Alar potent inhibitors lacking affinity for the co-factor PLP.     

Unconstrained enhanced sampling for free energy calculations of biomolecules: a review

Abstract

Free energy calculations are central to understanding the structure, dynamics and function of biomolecules. Yet insufficient sampling of biomolecular configurations is often regarded as one of the main sources of error. Many enhanced sampling techniques have been developed to address this issue. Notably, enhanced sampling methods based on biasing collective variables (CVs), including the widely used umbrella sampling, adaptive biasing force and metadynamics, have been discussed in a recent excellent review (Abrams and Bussi, Entropy, 2014). Here, we aim to review enhanced sampling methods that do not require predefined system-dependent CVs for biomolecular simulations and as such do not suffer from the hidden energy barrier problem as encountered in the CV-biasing methods. These methods include, but are not limited to, replica exchange/parallel tempering, self-guided molecular/Langevin dynamics, essential energy space random walk and accelerated molecular dynamics. While it is overwhelming to describe all details of each method, we provide a summary of the methods along with the applications and offer our perspectives. We conclude with challenges and prospects of the unconstrained enhanced sampling methods for accurate biomolecular free energy calculations.     

Study of the insulin dimerization: binding free energy calculations and per-residue free energy decomposition

A calculation of the binding free energy for the dimerization of insulin has been performed using the molecular mechanics-generalized Born surface area approach. The calculated absolute binding free energy is -11.9 kcal/mol, in approximate agreement with the experimental value of -7.2 kcal/mol. The results show that the dimerization is mainly due to nonpolar interactions. The role of the hydrogen bonds between the 2 monomers appears to give the direction of the interactions. A per-atom decomposition of the binding free energy has been performed to identify the residues contributing most to the self association free energy. Residues B24-B26 are found to make the largest favorable contributions to the dimerization. Other residues situated at the interface between the 2 monomers were found to make favorable but smaller contributions to the dimerization: Tyr B16, Val B12, and Pro B28, and to an even lesser extent, Gly B23. The energy decomposition on a per-residue basis is in agreement with experimental alanine scanning data. The results obtained from a single trajectory (i.e., the dimer trajectory is also used for the monomer analysis) and 2 trajectories (i.e., separate trajectories are used for the monomer and dimer) are similar.     

Estimation of the free energy of hard-sphere crystals via a free-volume approach

The free energies of the face-centred (FCC) and base-centred cubic (BCC) hard-sphere (HS) crystals have been estimated via the free-volume approach. We present two free-volume equations for the FCC and BCC HS crystals, which are different from those predicted by Velasco et al. [Langmuir 14(19) (1998), 5652–5655], and the equations exhibit more accuracy than Velasco et al.'s equations. The limitation of using the free-volume approach was assessed by comparing with true free energies obtained by other well-known methods, and possible improvement and application are discussed.     

Application of molecular dynamics and free energy perturbation methods to metalloporphyrin-ligand systems II: CO and dioxygen binding to myoglobin.

The protein contribution to the relative binding affinity of the ligands CO and O2 toward myoglobin (Mb) has been simulated using free energy perturbation calculations. The tautomers of the His E7 residue are different for the oxymyoglobin (MbO2) and carboxymyoglobin (MbCO) systems. This was modeled by performing two-step calculations that mutate the ligand and mutate the His E7 tautomers in separate steps. Differences in hydrogen bonding to the O2 and CO ligands were incorporated into the model. The O2 complex was calculated to be 2-3 kcal/mol more stable than the corresponding CO complex when compared to the same difference in an isolated heme control. This value agrees well with the experimental value of 2.0 kcal/mol. In qualitative agreement with experiments, the Fe-C-O bond is found to be bent (theta = 159.8 degrees) with a small tilt (theta = 6.2 degrees). The contributions made by each of the 29 residues--within the 9.0-A radius of the iron atom--to the free energy difference are separated into van der Waals and electrostatic contributions; the latter contributions are dominant. Aside from the proximal histidine and the heme group, the residues having the largest difference in free energy in mutating MbO2-->MbCO are His E7, Phe CD1, Phe CD4, Val E11, and Thr E10.     

Flavonoids as potent allosteric inhibitors of protein tyrosine phosphatase 1B: molecular dynamics simulation and free energy calculation

Protein tyrosine phosphatase 1B (PTP1B) is a member of the PTP superfamily which is considered to be a negative regulator of insulin receptor (IR) signaling pathway. PTP1B is a promising drug target for the treatment of type 2 diabetes, obesity, and cancer. The existence of allosteric site in PTP1B has turned the researcher’s attention to an alternate strategy for inhibition of this enzyme. Herein, the molecular interactions between the allosteric site of PTP1B with three non-competitive flavonoids, (MOR), (MOK), and (DPO) have been investigated. Three ligands were docked into allosteric site of the enzyme. The resulting protein–ligand complexes were used for molecular dynamics studies. Principal component and free-energy landscape (FEL) as well as cluster analyses were used to investigate the conformational and dynamical properties of the protein and identify representative enzyme substrates bounded to the inhibitors. Per residue energy decomposition analysis attributed dissimilar affinities of three inhibitors to the several hydrogen bonds and non-bonded interactions. In conclusion, our results exhibited an inhibitory pattern of the ligands against PTP1B.     

3D-QSAR-aided design of potent c-Met inhibitors using molecular dynamics simulation and binding free energy calculation

Abstract

Mesenchymal-epithelial transition factor (c-Met) is a member of receptor tyrosine kinase. It involves in various cellular signaling pathways which includes proliferation, motility, migration, and invasion. Over-expression of c-Met has been reported in various cancers. Hence, it is an ideal therapeutic target for cancer. The main objective of the study is to identify crucial residues involved in the inhibition of c-Met kinase and to design a series of potent imidazo [4,5-b] pyrazine derivatives as c-Met inhibitors. Docking was used to identify important active site residues involved in the inhibition of c-Met kinase which was further validated by 100 ns of molecular dynamics simulation and free energy calculation using molecular mechanics generalized born surface area. Furthermore, binding energy decomposition identified that residues Tyr1230, Met1211, Asp1222, Tyr1159, Met1160, Val1092, Ala1108, and Leu1157 contributed favorably to the binding stability of compound 32. Receptor-guided Comparative Molecular Field Analysis (CoMFA) (q² = 0.751, NOC = 6, r² = 0.933) and Comparative Molecular Similarity Indices Analysis (COMSIA) (q² = 0.744, NOC = 6, r² = 0.950) models were generated based on the docked conformation of the most active compound 32. The robustness of these models was tested using various validation techniques and found to be predictive. The results of CoMFA and CoMSIA contour maps exposed the regions favorable to enhance the activity. Based on this information, 27 novel c-Met inhibitors were designed. These designed compounds exhibited potent activity than the most active compound of the existing dataset.

Communicated by Ramaswamy H. Sarma     

Improved protein free energy calculation by more accurate treatment of nonbonded energy: Application to chymotrypsin inhibitor 2, V57A

We developed a software package for improved free energy calculation, in which spherical solvent boundary potential, cell multipole method, and Nosé-Hoover equation are employed. The performance of the developed software package is demonstrated in the case of valine to alanine mutation of the 57th residue in chymotrypsin inhibitor 2. By using this package, we obtained results quantitatively comparable to experimental results. By the free energy component analysis, it is shown that leucine 51, arginine 65, arginine 67, and phenylalanine 69 residues contribute significantly to the total free energy shift, ΔΔG. Among them, contribution from the hydrophilic arginine 67 residue, which is in close contact with the mutation site, is the largest. Structure around the mutation site is largely changed by the mutation. The structure change is caused mainly by two effects, hydrophobic interaction and short-range interaction along the sequence. Effects of Nosé-Hoover algorithm and Kirkwood reaction field are also discussed. Proteins 30:388–400, 1998. © 1998 Wiley-Liss, Inc.     

An explorative study on Staphylococcus aureus MurE inhibitor: induced fit docking,binding free energy calculation,and molecular dynamics

Staphylococcus aureus MurE enzyme catalyzes the addition of l-lysine as third residue of the peptidoglycan peptide moiety. Due to the high substrate specificity and its ubiquitous nature among bacteria, MurE enzyme is considered as one of the potential target for the development of new therapeutic agents. In the present work, induced fit docking (IFD), binding free energy calculation, and molecular dynamics (MD) simulation were carried out to elucidate the inhibition potential of 2-thioxothiazolidin-4-one based inhibitor 1 against S. aureus MurE enzyme. The inhibitor 1 formed majority of hydrogen bonds with the central domain residues Asn151, Thr152, Ser180, Arg187, and Lys219. Binding free-energy calculation by MM-GBSA approach showed that van der Waals (ΔG_vdW, ?57.30?kcal/mol) and electrostatic solvation (ΔG_solv, ?36.86?kcal/mol) energy terms are major contributors for the inhibitor binding. Further, 30-ns MD simulation was performed to validate the stability of ligand–protein complex and also to get structural insight into mode of binding. Based on the IFD and MD simulation results, we designed four new compounds D1–D4 with promising binding affinity for the S. aureus MurE enzyme. The designed compounds were subjected to the extra-precision docking and binding free energy was calculated for complexes. Further, a 30-ns MD simulation was performed for D1/4C13 complex.     

Surface free energy and interaction of Staphylococcus epidermidis with biomaterials

The adhesion of twenty nine Staphylococcus epidermidis strains to teflon, polyethylene, polycarbonate and bovine pericardium was studied in vitro and examined in relation to the surface free energies of both bacteria and biomaterials. All S. epidermidis strains had similar surface free energies, close to that of water, and adhered better to the materials with analogous surface free energies. There was a significant correlation (Kendall's Tau B = 1000) of biomaterial's surface free energy with the number of adhering bacteria. This correlation is inverse (Kendall's Tau B = -1000) when surface hydrophobicity is considered instead of surface free energy. This indicates that in Staphylococcus epidermidis adherence to biomaterials is inversely correlated to the surface hydrophobicity of the last, being so just the opposite of that occurring with other bacteria.