Computer simulations of ligand–protein binding with ensembles of protein conformations: A Monte Carlo study of HIV‐1 protease binding energy landscapes

We present the results of molecular docking simulations with HIV‐1 protease for the sb203386 and skf107457 inhibitors by Monte Carlo simulated annealing. A simplified piecewise linear energy function, the standard AMBER force field, and the AMBER force field with solvation and a soft‐core smoothing component are employed in simulations with a single‐protein conformation to determine the relationship between docking simulations with a simple energy function and more realistic force fields. The temperature‐dependent binding free energy profiles of the inhibitors interacting with a single protein conformation provide a detailed picture of relative thermodynamic stability and a distribution of ligand binding modes in agreement with experimental crystallographic data. Using the simplified piecewise linear energy function, we also performed Monte Carlo docking simulations with an ensemble of protein conformations employing preferential biased sampling of low‐energy protein conformations, and the results are analyzed in connection with the free energy profiles. ©1999 John Wiley & Sons, Inc. Int J Quant Chem 72: 73–84, 1999     

A conformational analysis method for understanding the energy landscapes of clusters.

A newly developed unbiased structural optimization method, named dynamic lattice searching (DLS), is proposed as an approach for conformational analysis of atomic/molecular clusters and used in understanding the energy landscape of large clusters. The structures of clusters are described in terms of the number of basic tetrahedron (BT) units they contain. We found that the hit numbers of different structural motifs in DLS runs is proportional to the number of BTs. A parameter T(max) is defined to limit the maximal number of atoms moved in a structural transition. Results show that T(max) is a key parameter for modulating the efficiency of the DLS method and has a great influence on the hit number of different motifs in DLS runs. Finally, the effect of potential range on the conformational distribution of the (Morse)(98) cluster is also discussed with different potential-range parameters.     

New optimization method for conformational energy calculations on polypeptides: Conformational space annealing

A new optimization method is presented to search for the global minimum-energy conformations of polypeptides. The method combines essential aspects of the build-up procedure and the genetic algorithm, and it introduces the important concept of “conformational space annealing.” Instead of considering a single conformation, attention is focused on a population of conformations while new conformations are obtained by modifying a “seed conformation.” The annealing is carried out by introducing a distance cutoff, D_cut, which is defined in the conformational space; D_cut effectively divides the whole conformational space of local minima into subdivisions. The value of D_cut is set to a large number at the beginning of the algorithm to cover the whole conformational space, and annealing is achieved by slowly reducing it. Many distinct local minima designed to be distributed as far apart as possible in conformational space are investigated simultaneously. Therefore, the new method finds not only the global minimum-energy conformation but also many other distinct local minima as by-products. The method is tested on Met-enkephalin, a 24-dihedral angle problem. For all 100 independent runs, the accepted global minimum-energy conformation was obtained after about 2600 minimizations on average. © 1997 John Wiley & Sons, Inc. J Comput Chem 18: 1222–1232     

Computing ensembles of transitions from stable states: Dynamic importance sampling

There is an increasing dataset of solved biomolecular structures in more than one conformation and increasing evidence that large-scale conformational change is critical for biomolecular function. In this article, we present our implementation of a dynamic importance sampling (DIMS) algorithm that is directed toward improving our understanding of important intermediate states between experimentally defined starting and ending points. This complements traditional molecular dynamics methods where most of the sampling time is spent in the stable free energy wells defined by these initial and final points. As such, the algorithm creates a candidate set of transitions that provide insights for the much slower and probably most important, functionally relevant degrees of freedom. The method is implemented in the program CHARMM and is tested on six systems of growing size and complexity. These systems, the folding of Protein A and of Protein G, the conformational changes in the calcium sensor S100A6, the glucose-galactose-binding protein, maltodextrin, and lactoferrin, are also compared against other approaches that have been suggested in the literature. The results suggest good sampling on a diverse set of intermediates for all six systems with an ability to control the bias and thus to sample distributions of trajectories for the analysis of intermediate states.     

Calculation of free energy landscapes: A histogram reweighted metadynamics approach

We present an efficient method for the calculation of free energy landscapes. Our approach involves a history‐dependent bias potential, which is evaluated on a grid. The corresponding free energy landscape is constructed via a histogram reweighting procedure a posteriori. Because of the presence of the bias potential, it can be also used to accelerate rare events. In addition, the calculated free energy landscape is not restricted to the actual choice of collective variables and can in principle be extended to auxiliary variables of interest without further numerical effort. The applicability is shown for several examples. We present numerical results for the alanine dipeptide and the Met‐Enkephalin in explicit solution to illustrate our approach. Furthermore, we derive an empirical formula that allows the prediction of the computational cost for the ordinary metadynamics variant in comparison with our approach, which is validated by a dimensionless representation. © 2011 Wiley Periodicals, Inc. J Comput Chem, 2011     

Basic ingredients of free energy calculations: A review

Methods to compute free energy differences between different states of a molecular system are reviewed with the aim of identifying their basic ingredients and their utility when applied in practice to biomolecular systems. A free energy calculation is comprised of three basic components: (i) a suitable model or Hamiltonian, (ii) a sampling protocol with which one can generate a representative ensemble of molecular configurations, and (iii) an estimator of the free energy difference itself. Alternative sampling protocols can be distinguished according to whether one or more states are to be sampled. In cases where only a single state is considered, six alternative techniques could be distinguished: (i) changing the dynamics, (ii) deforming the energy surface, (iii) extending the dimensionality, (iv) perturbing the forces, (v) reducing the number of degrees of freedom, and (vi) multi‐copy approaches. In cases where multiple states are to be sampled, the three primary techniques are staging, importance sampling, and adiabatic decoupling. Estimators of the free energy can be classified as global methods that either count the number of times a given state is sampled or use energy differences. Or, they can be classified as local methods that either make use of the force or are based on transition probabilities. Finally, this overview of the available techniques and how they can be best used in a practical context is aimed at helping the reader choose the most appropriate combination of approaches for the biomolecular system, Hamiltonian and free energy difference of interest. © 2009 Wiley Periodicals, Inc. J Comput Chem, 2010     

Conformational sensitivity of polyether macrocycles to electrostatic potential: Partial atomic charges,molecular mechanics,and conformational prediction

Molecular recognition (whether by enzymes, the immune system, or chelating ligands) depends critically on molecular conformation. Molecular mechanics predicts energetically favorable molecular conformations by locating low energy conformations using an empirical fit of molecular potential energy as a function of internal coordinates. Molecular mechanics analysis of 18-crown-6 demonstrates that the nonbonded term (primarily the electrostatic part) is the largest contributor to the conformational energy. Nevertheless, common methods of treating the electrostatic interaction for 18-crown-6 yield inconsistent values for conformational energies partly because partial charges assigned to each atom can change with conformation due to through-space inductive effects which are not considered in most molecular mechanics programs. Similar findings from several other groups are reviewed to support our conclusions. We argue for care and caution in predicting conformational preferences of molecules with two or more highly polar atoms. We also discuss the desirability of using an empirical method of partial charge determination such as the charge equilibration algorithm of Rappé and Goddard (or a suitable generalization which includes polarization) as a method of including these effects in molecular mechanics and molecular dynamics calculations.     

Systematic and statistical error in histogram-based free energy calculations

A common technique for the numerical calculation of free energies involves estimation of the probability density along a given coordinate from a set of configurations generated via simulation. The process requires discretization of one or more reaction coordinates to generate a histogram from which the continuous probability density is inferred. We show that the finite size of the intervals used to construct the histogram leads to quantifiable systematic error. The width of these intervals also determines the statistical error in the free energy, and the choice of the appropriate interval is therefore driven by the need to balance the two sources of error. We present a method for the construction of the optimal histogram for a given system, and show that the use of this technique requires little additional computational expense. We demonstrate the efficacy of the technique for a model system, and discuss how the principles governing the choice of discretization interval could be used to improve extended sampling techniques.     

键能的分子轨道理论研究 1: 理论公式

从LCAO-MO出发, 给出了一个计算键能的近似方法, 即EAB(i)-∑∑CaiSabCbiεi为第i个占据分子轨道(MO)中的一对电子对A-B键键能的贡献。对所有分子轨道求和即为该键的键能: EAB=∑EAB(i)。按该方法, 不仅可以计算各种不同分子中每两个相键连原子间的键能, 还可以从MO及AO角度分析每一具体键, 如σ, π, δ键的键能以及各AO对键能的贡献。该方法虽有别于求键焓和平衡离解能De, 但计算结果和De的实验值甚相符合。通过对键能的分析研究, 能较好地揭示原子间的相互作用关系及化学键的强弱, 从而可进一步探讨化学反应活性, 反应速率等化学性质。     

Ermod: Fast and versatile computation software for solvation free energy with approximate theory of solutions

ERmod is a software package to efficiently and approximately compute the solvation free energy using the method of energy representation. Molecular simulation is to be conducted at two condensed‐phase systems of the solution of interest and the reference solvent with test‐particle insertion of the solute. The subprogram ermod in ERmod then provides a set of energy distribution functions from the simulation trajectories, and another subprogram slvfe determines the solvation free energy from the distribution functions through an approximate functional. This article describes the design and implementation of ERmod, and illustrates its performance in solvent water for two organic solutes and two protein solutes. Actually, the free‐energy computation with ERmod is not restricted to the solvation in homogeneous medium such as fluid and polymer and can treat the binding into weakly ordered system with nano‐inhomogeneity such as micelle and lipid membrane. ERmod is available on web at http://sourceforge.net/projects/ermod . © 2014 Wiley Periodicals, Inc.     

Conformational Analysis of 2,2′-bifuran: Correlated High-level Ab initio and DFT Results

The torsional potential for inter-ring rotation in 2,2′-bifuran has been systematically tackled using highly accurate ab initio calculations as well as cost-effective DFT methods. The successful convergence of the ab initio results allowed to confirm the presence of a shallow gauche minimum in the torsional potential curve. The standard DFT methods failed to capture such a tiny energy barrier but, interestingly, the results could be remarkably improved by a mixture of wavefunction and DFT energies in a multi-coefficient fashion; thus, accurate DFT-based and ab initio reference data also become available. Since the experimental evaluation of torsional potentials faces quantitative problems, the outcome of high-level theoretical calculations is expected to be reliably used in further investigation on structure and conformational distribution of this system.     

Hybridizing rapidly exploring random trees and basin hopping yields an improved exploration of energy landscapes

The number of local minima of the potential energy landscape (PEL) of molecular systems generally grows exponentially with the number of degrees of freedom, so that a crucial property of PEL exploration algorithms is their ability to identify local minima, which are low lying and diverse. In this work, we present a new exploration algorithm, retaining the ability of basin hopping (BH) to identify local minima, and that of transition based rapidly exploring random trees (T‐RRT) to foster the exploration of yet unexplored regions. This ability is obtained by interleaving calls to the extension procedures of BH and T‐RRT, and we show tuning the balance between these two types of calls allows the algorithm to focus on low lying regions. Computational efficiency is obtained using state‐of‐the art data structures, in particular for searching approximate nearest neighbors in metric spaces. We present results for the BLN69, a protein model whose conformational space has dimension 207 and whose PEL has been studied exhaustively. On this system, we show that the propensity of our algorithm to explore low lying regions of the landscape significantly outperforms those of BH and T‐RRT. © 2015 Wiley Periodicals, Inc.     

键能的分子轨道理论研究2: 键能与Mulliken布居对键强度的 判断的比较

用键能E~A~B和Mulliken布居对化学键强度的判别进行了分析比较。结果表明,键能判据比Mulliken布居判据所得结论更符合实际情况。作为衡量原子间化学键强度的尺度,不仅应考虑原子轨道间的布居因素,还应考虑分子轨道(或原子轨道)的能量因素。     

Atropisomerism in Diarylamines: Structural Requirements and Mechanisms of Conformational Interconversion

In common with other hindered structures containing two aromatic rings linked by a short tether, diarylamines may exhibit atropisomerism (chirality due to restricted rotation). Previous examples have principally been tertiary amines, especially those with cyclic scaffolds. Little is known of the structural requirement for atropisomerism in structurally simpler secondary and acyclic diarylamines. In this paper we describe a systematic study of a series of acyclic secondary diarylamines, and we quantify the degree of steric hindrance in the ortho positions that is required for atropisomerism to result. Through a detailed experimental and computational analysis, the role of each ortho-substituent on the mechanism and rate of conformational interconversion is rationalised. We also present a simple predictive model for the design of configurationally stable secondary diarylamines.     

Conformational Analyses of Glycinal and Alaninal: A Computational Study

We present computational results from detailed gas-phase conformational analyses of the -substituted aldehydes, glycinal and alaninal. A synplanar conformer of glycinal and a synperiplanar conformer of alaninal in which the C=O and C–N bonds are in an eclipsing orientation are found to be lowest in energy; the two amino hydrogen atoms in these conformers are both directed over the C–C bond, i.e., in a compact arrangement. For the Group VA analogs, H₂P–CH₂–CHO and H₂P–CH(CH₃)–CHO, skew conformers in which the C–H and C–Me groups, respectively, are in an eclipsing orientation with the C=O bond are found to be lower in energy than the syn(peri)planar conformers. The results of various self-consistent reaction-field calculations suggest that the lowest-energy conformer of glycinal in 1, 2-dichloroethane is still synperiplanar, although the orientation of the amino hydrogen atoms may be different from that in the gas phase. Similar reaction-field calculations for alaninal raise the possibility that in this solvent a skew conformer, in which the C–H bond is nearly eclipsing the C=O bond, is energetically competitive with synperiplanar conformers.     

Conformational analysis of methylphenidate: comparison of molecular orbital and molecular mechanics methods

Summary Methylphenidate (MP) binds to the cocaine binding site on the dopamine transporter and inhibits reuptake of dopamine, but does not appear to have the same abuse potential as cocaine. This study, part of a comprehensive effort to identify a drug treatment for cocaine abuse, investigates the effect of choice of calculation technique and of solvent model on the conformational potential energy surface (PES) of MP and a rigid methylphenidate (RMP) analogue which exhibits the same dopamine transporter binding affinity as MP. Conformational analysis was carried out by the AM1 and AM1/SM5.4 semiempirical molecular orbital methods, a molecular mechanics method (Tripos force field with the dielectric set equal to that of vacuum or water) and the HF/6-31G* molecular orbital method in vacuum phase. Although all three methods differ somewhat in the local details of the PES, the general trends are the same for neutral and protonated MP. In vacuum phase, protonation has a distinctive effect in decreasing the regions of space available to the local conformational minima. Solvent has little effect on the PES of the neutral molecule and tends to stabilize the protonated species. The random search (RS) conformational analysis technique using the Tripos force field was found to be capable of locating the minima found by the molecular orbital methods using systematic grid search. This suggests that the RS/Tripos force field/vacuum phase protocol is a reasonable choice for locating the local minima of MP. However, the Tripos force field gave significantly larger phenyl ring rotational barriers than the molecular orbital methods for MP and RMP. For both the neutral and protonated cases, all three methods found the phenyl ring rotational barriers for the RMP conformers/invertamers (denoted as cte, tte, and cta) to be: cte, tte> MP > cta. Solvation has negligible effect on the phenyl ring rotational barrier of RMP. The B3LYP/6-31G* density functional method was used to calculate the phenyl ring rotational barrier for neutral MP and gave results very similar to those of the HF/6-31G* method.     

Analytical potential energy functions and theoretical spectroscopic constants for MX/MX− (MGe,Sn, Pb; XO,S, Se,Te, Po) and LuA (AH,F) systems: Density functional theory calculations

The molecular spectroscopic constants for the chalcogenide complexes MX (M?Ge, Sn, Pb; X?O, S, Se, Te, Po) and their corresponding MX^? anions are presented with the LSDA/SDD, B1LYP/SDD, and B3LYP/SDD methods. Although many methods were attempted, only the most promising results are reported. We show that the best results are obtained by LSDA/SDD calculations, and thus this method is emphasized as an illustrative example of our methodology. The potential energy curves and physical property characterizations for X¹ ∑⁺ state of LuH and LuF are presented with a variety of density functional theory (DFT) methods. Comparisons with wave function‐based treatments (HF, MP2, CCSD, QCISD) are made in addition to experimental correlations. We show that the best results are obtained by the B3LYP/SDD method for LuH, and the MPW1PW91/SDD method for LuF. © 2006 Wiley Periodicals, Inc. Int J Quantum Chem, 2007     

Exploring potential energy surfaces for chemical reactions: an overview of some practical methods

Potential energy surfaces form a central concept in the application of electronic structure methods to the study of molecular structures, properties, and reactivities. Recent advances in tools for exploring potential energy surfaces are surveyed. Methods for geometry optimization of equilibrium structures, searching for transition states, following reaction paths and ab initio molecular dynamics are discussed. For geometry optimization, topics include methods for large molecules, QM/MM calculations, and simultaneous optimization of the wave function and the geometry. Path optimization methods and dynamics based techniques for transition state searching and reaction path following are outlined. Developments in the calculation of ab initio classical trajectories in the Born-Oppenheimer and Car-Parrinello approaches are described.     

Energy landscapes of nucleophilic substitution reactions: a comparison of density functional theory and coupled cluster methods

We have carried out a detailed evaluation of the performance of all classes of density functional theory (DFT) for describing the potential energy surface (PES) of a wide range of nucleophilic substitution (SN2) reactions involving, amongst others, nucleophilic attack at carbon, nitrogen, silicon, and sulfur. In particular, we investigate the ability of the local density approximation (LDA), generalized gradient approximation (GGA), meta-GGA as well as hybrid DFT to reproduce high-level coupled cluster (CCSD(T)) benchmarks that are close to the basis set limit. The most accurate GGA, meta-GGA, and hybrid functionals yield mean absolute deviations of about 2 kcal/mol relative to the coupled cluster data, for reactant complexation, central barriers, overall barriers as well as reaction energies. For the three nonlocal DFT classes, the best functionals are found to be OPBE (GGA), OLAP3 (meta-GGA), and mPBE0KCIS (hybrid DFT). The popular B3LYP functional is not bad but performs significantly worse than the best GGA functionals. Furthermore, we have compared the geometries from several density functionals with the reference CCSD(T) data. The same GGA functionals that perform best for the energies (OPBE, OLYP), also perform best for the geometries with average absolute deviations in bond lengths of 0.06 A and 0.6 degrees, even better than the best meta-GGA and hybrid functionals. In view of the reduced computational effort of GGAs with respect to meta-GGAs and hybrid functionals, let alone coupled cluster, we recommend the use of accurate GGAs such as OPBE or OLYP for the study of SN2 reactions.