Parallelized-over-parts computation of absolute binding free energy with docking and molecular dynamics

We present a technique for biomolecular free energy calculations that exploits highly parallelized sampling to significantly reduce the time to results. The technique combines free energies for multiple, nonoverlapping configurational macrostates and is naturally suited to distributed computing. We describe a methodology that uses this technique with docking, molecular dynamics, and free energy perturbation to compute absolute free energies of binding quickly compared to previous methods. The method does not require a priori knowledge of the binding pose as long as the docking technique used can generate reasonable binding modes. We demonstrate the method on the protein FKBP12 and eight of its inhibitors.     

Extra precision docking,free energy calculation and molecular dynamics studies on glutamic acid derivatives as MurD inhibitors

The binding modes of well known MurD inhibitors have been studied using molecular docking and molecular dynamics (MD) simulations. The docking results of inhibitors 1-30 revealed similar mode of interaction with Escherichia coli-MurD. Further, residues Thr36, Arg37, His183, Lys319, Lys348, Thr321, Ser415 and Phe422 are found to be important for inhibitors and E. coli-MurD interactions. Our docking procedure precisely predicted crystallographic bound inhibitor 7 as evident from root mean square deviation (0.96 Å). In addition inhibitors 2 and 3 have been successfully cross-docked within the MurD active site, which was pre-organized for the inhibitor 7. Induced fit best docked poses of 2, 3, 7 and 15/2Y1O complexes were subjected to 10 ns MD simulations to determine the stability of the predicted binding conformations. Induce fit derived docked complexes were found to be in a state of near equilibrium as evident by the low root mean square deviations between the starting complex structure and the energy minimized final average MD complex structures. The results of molecular docking and MD simulations described in this study will be useful for the development of new MurD inhibitors with high potency.     

3D-QSAR,docking, molecular dynamics simulation and free energy calculation studies of some pyrimidine derivatives as novel JAK3 inhibitors

Janus kinase 3 (JAK3) is a promising drug target for the treatment of inflammatory diseases, autoimmune disorders, organ transplant rejection and various cancers. In the present study, 3D-QSAR, docking, MD simulation and MM/PBSA studies were performed on a series of pyrimidine-based JAK3 inhibitors. A reliable COMSIA (q² = 0.717 and r² = 0.986) model was developed and validated using external validation test set, bootstrapping, progressive scrambling and r_m² metrics analyses. Structural requirements identified through contour maps of the model were strategically utilized to computationally design 170 novel JAK3 inhibitors with improved potency. Docking studies were performed on the selected data set and newly designed compounds to show their binding mode and to identify important interacting residues inside the active site of JAK3. In addition, docking results of the selected designed compounds inside the active sites of JAK1, JAK2 and TYK2 indicated their JAK3 selectivity. MD simulation (100 ns) on the docked complex of compound 28 (one of highly active compounds of the data set) assisted in the further exploration of the binding interactions. Some crucial residues like Lys830 (glycine-rich loop), Val836, Ala853, Leu905 (hinge region), Cys909, Asn954, Leu956 and Ala966 were identified. Hydrogen bond interactions with hinge residue Leu905 were critical for the binding of JAK3 inhibitors. Additionally, MM/PBSA calculation provided the binding free energy of the compound 28. Newly designed molecules showed promising results in the preliminary in silico ADMET evaluations. Outcomes of the study can further be exploited to develop potent JAK3 inhibitors.     

Direct calculation of the crystal-melt interfacial free energy via molecular dynamics computer simulation

We review our recent work on the direct calculation of the interfacial free energy, gamma, of the crystal-melt interface via molecular dynamics computer simulation for a number of model systems. The value of gamma as a function of crystal orientation is determined using a thermodynamic integration technique employing moving cleaving walls [Phys. Rev. Lett. 2000, 85, 4751]. The calculation is sufficiently precise to resolve the small anisotropy in gamma, which is crucial in determining the kinetics and morphology of dendritic growth. We report values of gamma for the hard-sphere and Lennard-Jones systems, as well as recent results on the series of inverse-power potentials. For the inverse sixth-, seventh-, and eighth-power systems, we determine gamma for both fcc and bcc crystal structures. For these systems, the bcc-melt gamma is lower than that for fcc crystals by about 25%, consistent with recent experiments and computer simulations on fcc-forming systems that show preferential formation of bcc nuclei in the initial stages of crystallization. In addition, we show that our results give a molecular interpretation of Turnbull's rule, which is the empirical relationship between gamma and the enthalpy of fusion.     

Phase diagrams of alkali halides using two interaction models: a molecular dynamics and free energy study

Phase diagrams for potassium and sodium chlorides are determined by molecular dynamics and free energy calculations. Two rigid-ion interaction models, namely, the Born-Mayer-Huggins (BMH) and Michielsen-Woerlee-Graaf (MWG) effective pair potentials, have been used. The critical and triple point properties are discussed and compared with available experimental and simulation data. The MWG model reproduces the experimental liquid-gas equilibria better than the BMH model, being the accordance very good in the lowest temperature region of the coexistent liquids, particularly for NaCl. However, both models underestimate the critical temperatures of KCl and NaCl. Relatively to the solid-gas equilibria, the models do not reproduce well the experimental data. As for the solid-liquid coexistences either the BMH or the MWG models appear unrealistic.     

Medicinal plant phytochemicals and their inhibitory activities against pancreatic lipase: molecular docking combined with molecular dynamics simulation approach

Obesity is the worst health risk worldwide, which is linked to a number of diseases. Pancreatic lipase is considered as an affective cause of obesity and can be a major target for controlling the obesity. The present study was designed to find out best phytochemicals against pancreatic lipase through molecular docking combined with molecular dynamics (MD) simulation. For this purpose, a total of 3770 phytochemicals were docked against pancreatic lipase and ranked them on the basis of binding affinity. Finally, 10 molecules (Kushenol K, Rosmarinic acid, Reserpic acid, Munjistin, Leachianone G, Cephamycin C, Arctigenin, 3-O-acetylpadmatin, Geniposide and Obtusin) were selected that showed strong bonding with the pancreatic lipase. MD simulations were performed on top five compounds using AMBER16. The simulated complexes revealed stability and ligands remained inside the binding pocket. This study concluded that these finalised molecules can be used as drug candidate to control obesity.     

Sampling enhancement for the quantum mechanical potential based molecular dynamics simulations: a general algorithm and its extension for free energy calculation on rugged energy surface

An approach is developed in the replica exchange framework to enhance conformational sampling for the quantum mechanical (QM) potential based molecular dynamics simulations. Importantly, with our enhanced sampling treatment, a decent convergence for electronic structure self-consistent-field calculation is robustly guaranteed, which is made possible in our replica exchange design by avoiding direct structure exchanges between the QM-related replicas and the activated (scaled by low scaling parameters or treated with high "effective temperatures") molecular mechanical (MM) replicas. Although the present approach represents one of the early efforts in the enhanced sampling developments specifically for quantum mechanical potentials, the QM-based simulations treated with the present technique can possess the similar sampling efficiency to the MM based simulations treated with the Hamiltonian replica exchange method (HREM). In the present paper, by combining this sampling method with one of our recent developments (the dual-topology alchemical HREM approach), we also introduce a method for the sampling enhanced QM-based free energy calculations.     

Computer graphics in real-time docking with energy calculation and minimization

We describe a real-time docking method using molecular graphics and high-speed calculation of the energy of interaction between a drug and a receptor. A three-dimensional tabulation of the potential is computed prior to the docking experiment, and the total interaction energy is calculated and updated on the display screen in real-time to provide immediate visual feedback to the user as the drug molecule is moved inside the receptor pocket. The “Simplex” method is then used to minimize the energy of interactions after each docking. Using this real-time method, it is now possible to examine rapidly the interactions of a large number of drugs and their analogues with receptor molecules. As an application, the interaction of thyroid hormone and its analogues with prealbumin is considered. The final interaction energies of this very rapid method compare well with those calculated by more orthodox means, while also providing visual feedback on both molecular geometry and energy.     

Virtual screening of B-Raf kinase inhibitors: A combination of pharmacophore modelling,molecular docking, 3D-QSAR model and binding free energy calculation studies

B-Raf kinase has been identified as an important target in recent cancer treatment. In order to discover structurally diverse and novel B-Raf inhibitors (BRIs), a virtual screening of BRIs against ZINC database was performed by using a combination of pharmacophore modelling, molecular docking, 3D-QSAR model and binding free energy (ΔG_bind) calculation studies in this work. After the virtual screening, six promising hit compounds were obtained, which were then tested for inhibitory activities of A375 cell lines. In the result, five hit compounds show good biological activities (IC₅₀ < 50 μM). The present method of virtual screening can be applied to find structurally diverse inhibitors, and the obtained five structurally diverse compounds are expected to develop novel BRIs.     

A molecular dynamics study of the lateral free energy profile of a pair of cholesterol molecules as a function of their distance in phospholipid bilayers

Free energy profile of a pair of cholesterol molecules in a leaflet of 1-palmitoyl-2-oleoyl-phosphatidylcholine (POPC) bilayers in the liquid-crystalline phase has been calculated as a function of their lateral distance using a combination of NPT-constant atomistic molecular dynamics calculations (P = 1 atm and T = 310.15 K) and the thermodynamic integration method. The calculated free energy clearly shows that the two cholesterol molecules form a dimer separated by a distance of 1.0-1.5 nm in POPC bilayers. Well depth of the free energy profile is about 3.5 kJ/mol, which is comparable to the thermal energy k(B)T at 310.15 K. This indicates that the aggregation of cholesterol molecules in the bilayers depends on the temperature as well as the concentration of the system. The free energy function obtained here may be used as a reference when coarse grained potential model is investigated for this two-component system. Local structure of POPC molecules around two cholesterol molecules has also been investigated.     

Discovery of potential inhibitor against human acetylcholinesterase: a molecular docking and molecular dynamics investigation

Alzheimer’s disease (AD) is a progressive neurodegenerative disease of central nervous system among elderly people. Human acetylcholinesterase (hAChE), an important enzyme in neuronal signaling, is responsible for the degradation of acetylcholine which in turn prevents the post synaptic signal transmissions. hAChE has been an attractive target of drug discovery for the search of therapeutics against AD. In the recent past hAChE has become hot target for the investigation of new potential therapeutics. We performed virtual screening of entire database against hAChE. Further, the extra precision molecular docking was carried out to refine the docking results and the best complex was passed for molecular dynamics simulations in order of understanding the hAChE dynamics and its behavior in complex with the ligand which corroborate the outcomes of virtual screening. This also provides binding free energy data that establishes the ligands efficiency for inhibiting hAChE. The computational findings discussed in this paper provide initial information of inhibitory effects of ligand, (drugbank entry DB00983), over hAChE.     

Computation of binding free energy with molecular dynamics and grand canonical Monte Carlo simulations

The binding of a ligand to a receptor is often associated with the displacement of a number of bound water molecules. When the binding site is exposed to the bulk region, this process may be sampled adequately by standard unbiased molecular dynamics trajectories. However, when the binding site is deeply buried and the exchange of water molecules with the bulk region may be difficult to sample, the convergence and accuracy in free energy perturbation (FEP) calculations can be severely compromised. These problems are further compounded when a reduced system including only the region surrounding the binding site is simulated. To address these issues, we couple molecular dynamics (MD) with grand canonical Monte Carlo (GCMC) simulations to allow the number of water to fluctuate during an alchemical FEP calculation. The atoms in a spherical inner region around the binding pocket are treated explicitly while the influence of the outer region is approximated using the generalized solvent boundary potential (GSBP). At each step during thermodynamic integration, the number of water in the inner region is equilibrated with GCMC and energy data generated with MD is collected. Free energy calculations on camphor binding to a deeply buried pocket in cytochrome P450cam, which causes about seven water molecules to be expelled, are used to test the method. It concluded that solvation free energy calculations with the GCMC/MD method can greatly improve the accuracy of the computed binding free energy compared to simulations with fixed number of water.     

Minimum free energy pathways and free energy profiles for conformational transitions based on atomistic molecular dynamics simulations

An efficient method for the calculation of minimum free energy pathways and free energy profiles for conformational transitions is presented. Short restricted perturbation-targeted molecular dynamics trajectories are used to generate an approximate free energy surface. Approximate reaction pathways for the conformational change are constructed from one-dimensional line segments on this surface using a Monte Carlo optimization. Accurate free energy profiles are then determined along the pathways by means of one-dimensional adaptive umbrella sampling simulations. The method is illustrated by its application to the alanine "dipeptide." Due to the low computational cost and memory demands, the method is expected to be useful for the treatment of large biomolecular systems.     

The docking of chiral epoxides on the Whelk-O1 stationary phase: a molecular dynamics study

The docking of analytes on the Whelk-O1 chiral stationary phase is explored for two chiral epoxides in a hexane solvent. Density functional theory calculations are employed to develop flexible models for R/S-styrene oxide (phenyl oxirane) and (R,R/S,S)-stilbene oxide (2,3-diphenyl oxirane). Molecular dynamics simulations of the racemates in the presence of the Whelk-O1 chiral stationary phase reveal the distribution of the enantiomers at the interface. The importance of hydrogen bonding and ring-ring interactions is explored along with an examination of the major docking arrangements. The interactions between the Whelk-O1 molecules and the chiral epoxide enantiomers are quite distinct and consistent with the experimental elution orders [S.E. Schaus, B.D. Brandes, J.F. Larrow, M. Tokunage, K.B. Hansen, A.E. Gould, M.E. Furrow, E.N. Jacobsen, J. Am. Chem. Soc. 124 (2001) 1307] and separation factors [W.H. Pirkle, C.J. Welch, Tetrahedron: Asymm. 5 (1994) 777]. The impact of a polar solvent modifier is examined for R/S-styrene oxide where selectivity in 80:20 n-hexane:2-propanol is assessed.     

Variational calculation of quantum mechanical/molecular mechanical free energy with electronic polarization of solvent

Quantum mechanical/molecular mechanical (QM/MM) free energy calculation presents a significant challenge due to an excessive number of QM calculations. A useful approach for reducing the computational cost is that based on the mean field approximation to the QM subsystem. Here, we describe such a mean-field QM/MM theory for electronically polarizable systems by starting from the Hartree product ansatz for the total system and invoking a variational principle of free energy. The MM part is then recast to a classical polarizable model by introducing the charge response kernel. Numerical test shows that the potential of mean force (PMF) thus obtained agrees quantitatively with that obtained from a direct QM/MM calculation, indicating the utility of self-consistent mean-field approximation. Next, we apply the obtained method to prototypical reactions in several qualitatively different solvents and make a systematic comparison of polarization effects. The results show that in aqueous solution the PMF does not depend very much on the water models employed, while in nonaqueous solutions the PMF is significantly affected by explicit polarization. For example, the free energy barrier for a phosphoryl dissociation reaction in acetone and cyclohexane is found to increase by more than 10 kcal/mol when switching the solvent model from an empirical to explicitly polarizable one. The reason for this is discussed based on the parametrization of empirical nonpolarizable models.     

Some practical aspects of free energy calculations from molecular dynamics simulation

In this work, we address two critical aspects of calculation of the free energy differences in molecular systems from molecular simulations. The first aspect involves checking whether the calculated free energy difference depends significantly on the extent of perturbation used for accomplishment of a given transformation. The second aspect of interest is to verify if the sampling errors in calculating the free energy differences between the wild-type molecule and a mutated one in its free state and in a complex are similar, or not, for a finite-length dynamic simulation. The reliability of the free energy estimates obtained from molecular simulations using thermodynamic cycles depends in part on this fact. For investigating these aspects, we use a self-transformation scheme in which a transformation of a part of a molecular system into itself is considered. We perform MD simulations of DNA fragments in which a part of a specific base is subjected to such a self-transformation. Results indicate that the estimated free energy differences do not depend significantly on the extent of perturbation used to achieve the transformation. Interestingly, the variation in the cumulative free energy difference, ΔA, with the coupling parameter, λ, depends significantly on the extent of perturbation. We examine the physical basis of the observed nature of the variation of the accumulated free energy difference, ΔA, against the λ value in the case of a self-transformation. In a thermodynamic cycle, the sampling errors due to the finite-length simulation for the molecular system are found to be similar to each other for the two perturbations (free and in a complex) justifying the use of such approach in calculating ΔΔA in molecular complexes. ©1999 John Wiley & Sons, Inc. J Comput Chem 20: 877–885, 1999     

Efficient calculation of many stacking and pairing free energies in DNA from a few molecular dynamics simulations

Through the use of the one-step perturbation approach, 130 free energies of base stacking and 1024 free energies of base pairing in DNA have been calculated from only five simulations of a nonphysical reference state. From analysis of a diverse set of 23 natural and unnatural bases, it appears that stacking free energies and stacking conformations play an important role in pairing of DNA nucleotides. On the one hand, favourable pairing free energies were found for bases that do not have the possibility to form canonical hydrogen bonds, while on the other hand, good hydrogen-bonding possibilities do not guarantee a favourable pairing free energy if the stacking of the bases dictates an unfavourable conformation. In this application, the one-step perturbation approach yields a wealth of both energetic and structural information at minimal computational cost.     

Computational studies on pseudorotaxanes by molecular dynamics and free energy perturbation simulations

A computational scheme that comprises the utilization of the AMBER force field with RESP charges and an explicit solvent model for acetonitrile proved to be useful for studying the structures and energetics of pseudorotaxanes of benzidine and 4,4'-biphenol with cyclobis(paraquat-p-phenylene). The scheme can be further utilized for modeling [2]rotaxanes.     

基于分子对接的thanatin与NDM-1分子动力学模拟研究

耐碳青霉烯类抗生素的超级细菌给人类健康带来了严重威胁,其所携带的金属 β-内酰胺酶编码基因是耐药性的主要来源。NDM-1作为其中传播最广、活性最强的 β-内酰胺酶,其抑制剂的研发刻不容缓。具有广谱作用的抗菌肽thanatin对NDM-1展现出了较好的抑制效果,但抑制机理并不清楚。本文使用HPEPDOCK与Rosetta FlexPepDock服务器,将thanatin与NDM-1进行了分子对接,并使用Desmond软件包对对接模型进行了分子动力学模拟。结果表明,thanatin与NDM-1活性中心的Zn2+ 并无直接相互作用,而作为竞争性抑制剂结合于NDM-1的活性口袋,阻止抗生素分子进入活性口袋与Zn2+ 结合,从而抑制NDM-1的水解活性。本文为研发有效的NDM临床抑制剂探索了可行的方法。