Stochastic molecular optimization using generalized simulated annealing

We propose a stochastic optimization technique based on a generalized simulated annealing (GSA) method for mapping minima points of molecular conformational energy surfaces. The energy maps are obtained by coupling a classical molecular force field (THOR package) with a GSA procedure. Unlike the usual molecular dynamics (MD) method, the method proposed in this study is force independent; that is, we obtain the optimized conformation without calculating the force, and only potential energy is involved. Therefore, we do not need to know the conformational energy gradient to arrive at equilibrium conformations. Its utility in molecular mechanics is illustrated by applying it to examples of simple molecules (H₂O and H₂O₃) and to polypeptides. The results obtained for H₂O and H₂O₃ using Tsallis thermostatistics suggest that the GSA approach is faster than the other two conventional methods (Boltzmann and Cauchy machines). The results for polypeptides show that pentalanine does not form a stable α-helix structure, probably because the number of hydrogen bonds is insufficient to maintain the helical array. On the contrary, the icoalanine molecule forms an α-helix structure. We obtain this structure simulating all Φ, Ψ pairs using only a few steps, as compared with conventional methods. © 1998 John Wiley & Sons, Inc. J Comput Chem 19: 647–657, 1998     

Efficiency of the multicanonical simulation method as applied to peptides of increasing size: the heptapeptide deltorphin

The advantage of the multicanonical (MUCA) simulation method of Berg and coworkers over the conventional Metropolis method is in its ability to move a system effectively across energy barriers thereby providing results for a wide range of temperatures. However, a MUCA simulation is based on weights (related to the density of states) that should be determined prior to a production run and their calculation is not straightforward. To overcome this difficulty a procedure has been developed by Berg that calculates the MUCA weights automatically. In a previous article (Ya?ar et al. J Comput Chem 2000, 14, 1251-1261) we extended this procedure to continuous systems and applied it successfully to the small pentapeptide Leu-enkephalin. To investigate the performance of the automated MUCA procedure for larger peptides, we apply it here to deltorphin, a linear heptapeptide with bulky side chains (H-Tyr(1)-D-Met(2)-Phe(3)-His(4)-Leu(5)-Met(6)-Asp(7)-NH(2)). As for Leu-enkephalin, deltorphin is modeled in vacuum by the potential energy function ECEPP. MUCA is found to perform well. A weak second peak is seen for the specific heat, which is given a special attention. By minimizing the energy of structures along the trajectory it is found that MUCA provides a good conformational coverage of the low energy region of the molecule. These latter results are compared with conformational coverage obtained by the Monte Carlo minimization method of Li and Scheraga.     

The flexible alignment of molecular structures using simulated annealing with weighted Lagrangian multipliers

A framework for superimposing small molecules is presented. The proposed method consists of a simple atom‐based, flexible alignment. The optimization procedure used in the alignment is based on a recently published variant of the simulated annealing whereby nonlinear constraints are accommodated using Lagrangian multipliers. It differs from other published superposition algorithms in that any number of nonlinear constraints can be readily imposed on the structural alignment directly through the objective function without assuming an a priori trade‐off between competing conditions. These can include equality and equality constraints on distances, angles, and energy states. Examples illustrating the use of the proposed approach are also provided. © 2010 Wiley Periodicals, Inc. J Comput Chem, 2011     

Protein structure prediction aided by geometrical and probabilistic constraints

Database-assisted ab initio protein structure prediction methods have exhibited considerable promise in the recent past, with several implementations being successful in community-wide experiments (CASP). We have employed combinatorial optimization techniques toward solving the protein structure prediction problem. A Monte Carlo minimization algorithm has been employed on a constrained search space to identify minimum energy configurations. The search space is constrained by using radius of gyration cutoffs, the loop backbone dihedral probability distributions, and various secondary structure packing conformations. Simulations have been carried out on several sequences and 1000 conformations have been initially generated. Of these, 50 best candidates have then been selected as probable conformations. The search for the optimum has been simplified by incorporating various geometrical constraints on secondary structural elements using distance restraint potential functions. The advantages of the reported methodology are its simplicity, and modifiability to include other geometric and probabilistic restraints.     

Global optimization using bad derivatives: Derivative-free method for molecular energy minimization

A general method designed to isolate the global minimum of a multidimensional objective function with multiple minima is presented. The algorithm exploits an integral “coarse-graining” transformation of the objective function, U, into a smoothed function with few minima. When the coarse-graining is defined over a cubic neighborhood of length scale ϵ, the exact gradient of the smoothed function, 𝒰_ϵ, is a simple three-point finite difference of U. When ϵ is very large, the gradient of 𝒰_ϵ appears to be a “bad derivative” of U. Because the gradient of 𝒰_ϵ is a simple function of U, minimization on the smoothed surface requires no explicit calculation or differentiation of 𝒰_ϵ. The minimization method is “derivative-free” and may be applied to optimization problems involving functions that are not smooth or differentiable. Generalization to functions in high-dimensional space is straightforward. In the context of molecular conformational optimization, the method may be used to minimize the potential energy or, preferably, to maximize the Boltzmann probability function. The algorithm is applied to conformational optimization of a model potential, Lennard–Jones atomic clusters, and a tetrapeptide. © 1998 John Wiley & Sons, Inc. J Comput Chem 19: 1445–1455, 1998     

Comparative analysis of the conformational profile of substance P using simulated annealing and molecular dynamics

The present study describes an extensive conformational search of substance P using two different computational methods. On the one hand, the peptide was studied using the iterative simulated annealing, and on the other, molecular dynamics simulations at 300 and 400 K. With the former method, the peptide was studied in vacuo with a dielectric constant of 80, whereas using the latter study the peptide was studied in a box of TIP3P water molecules. Analysis of the results obtained using both methodologies was carried out using an in-house methodology using a cluster analysis method based on information theory. Comparison of the two sampling methodologies and the different environment used in the calculations is also analyzed. Finally, the conformational motifs that are characteristic of substance P in a hydrophilic environment are presented and compared with the experimental results available in the literature.     

Exploring the similarities between potential smoothing and simulated annealing

Simulated annealing and potential function smoothing are two widely used approaches for global energy optimization of molecular systems. Potential smoothing as implemented in the diffusion equation method has been applied to study partitioning of the potential energy surface (PES) for N‐Acetyl‐Ala‐Ala‐N‐Methylamide (CDAP) and the clustering of conformations on deformed surfaces. A deformable version of the united‐atom OPLS force field is described, and used to locate all local minima and conformational transition states on the CDAP surface. It is shown that the smoothing process clusters conformations in a manner consistent with the inherent structure of the undeformed PES. Smoothing deforms the original surface in three ways: structural shifting of individual minima, merging of adjacent minima, and energy crossings between unrelated minima. A master equation approach and explicit molecular dynamics trajectories are used to uncover similar features in the equilibrium probability distribution of CDAP minima as a function of temperature. Qualitative and quantitative correlations between the simulated annealing and potential smoothing approaches to enhanced conformational sampling are established. © 2000 John Wiley & Sons, Inc. J Comput Chem 21: 531–552, 2000     

An efficient molecular docking using conformational space annealing

Molecular docking falls into the general category of global optimization problems because its main purpose is to find the most stable complex consisting of a receptor and its ligand. Conformational space annealing (CSA), a powerful global optimization method, is incorporated with the Tinker molecular modeling package to perform molecular docking simulations of six receptor-ligand complexes (3PTB, 1ULB, 2CPP, 1STP, 3CPA, and 1PPH) from the Protein Data Bank. In parallel, Monte Carlo with the minimization (MCM) method is also incorporated into the Tinker package for comparison. The energy function, consisting of electrostatic interactions, van der Waals interactions, and torsional energy terms, is calculated using the AMBER94 all-atom empirical force field. Rigid docking simulations for all six complexes and flexible docking simulations for three complexes (1STP, 3CPA, and 1PPH) are carried out using the CSA and the MCM methods. The simulation results show that the docking procedures using the CSA method generally find the most stable complexes as well as the native-like complexes more efficiently and accurately than those using the MCM, demonstrating that CSA is a promising search method for molecular docking problems.     

Efficiency of tabu‐search‐based conformational search algorithms

Efficient conformational search or sampling approaches play an integral role in molecular modeling, leading to a strong demand for even faster and more reliable conformer search algorithms. This article compares the efficiency of a molecular dynamics method, a simulated annealing method, and the basin hopping (BH) approach (which are widely used in this field) with a previously suggested tabu‐search‐based approach called gradient only tabu search (GOTS). The study emphasizes the success of the GOTS procedure and, more importantly, shows that an approach which combines BH and GOTS outperforms the single methods in efficiency and speed. We also show that ring structures built by a hydrogen bond are useful as starting points for conformational search investigations of peptides and organic ligands with biological activities, especially in structures that contain multiple rings. © 2011 Wiley Periodicals, Inc. J Comput Chem, 2011     

Multicanonical procedure for continuum peptide models

The multicanonical (Muca) Monte Carlo method enables simulating a system over a wide range of temperatures and thus has become an efficient tool for studying spin glasses, first‐order phase transitions, the helix–coil transition of polypeptides, and protein folding. However, implementation of the method requires calculating the multicanonical weights by an iterative procedure that is not straightforward and is a stumbling block for newcomers. A recursive procedure that takes into account the statistical errors of all previous iterations and thus enables an automatic calculation of the weights without the need for human intervention after each iteration has been proposed. This procedure, which has already been tested successfully for lattice systems, is extended here to continuum models of peptides and proteins. The method is examined in detail and tested for models of the pentapeptide Leu‐enkephalin (Tyr‐Gly‐Gly‐Phe‐Leu) described by the potential energy function ECEPP. Because of the great interest in the structural mapping of the low‐energy region of biomolecules, the energy of structures selected from the Muca trajectory is minimized. The extent of conformational coverage provided by the method is examined and found to be very satisfactory. © 2000 John Wiley & Sons, Inc. J Comput Chem 21: 1251–1261, 2000     

Conformational search of peptides and proteins: Monte Carlo minimization with an adaptive bias method applied to the heptapeptide deltorphin

The energy function of a protein consists of a tremendous number of minima. Locating the global energy minimum (GEM) structure, which corresponds approximately to the native structure, is a severe problem in global optimization. Recently we have proposed a conformational search technique based on the Monte Carlo minimization (MCM) method of Li and Scheraga, where trial dihedral angles are not selected at random within the range [-180 degrees,180 degrees ] (as with MCM) but with biased probabilities depending on the increased structure-energy correlations as the GEM is approached during the search. This method, called the Monte Carlo minimization with an adaptive bias (MCMAB), was applied initially to the pentapeptide Leu-enkephalin. Here we study its properties further by applying it to the larger peptide with bulky side chains, deltorphin (H-Tyr-D-Met-Phe-His-Leu-Met-Asp-NH(2)). We find that on average the number of energy minimizations required by MCMAB to locate the GEM for the first time is smaller by a factor of approximately three than the number required by MCM-in accord with results obtained for Leu-enkephalin.     

Conformational behavior of antineoplastic peptides Dolastatin 10 and Dolastatin 15 from Monte Carlo and molecular dynamics simulations

In this work, we report the results of a combined Monte Carlo (MC) and molecular dynamics (MD) approach applied to the conformational study of natural potent mitosis inhibitor Dolastatin 15. A GRID and CLogP analysis has been performed with the aim to provide some hints to the interpretation of the different behavior of Dolastatin 15 and Dolastatin 10, which has been examined in a former study. © 2006 Wiley Periodicals, Inc. Int J Quantum Chem, 2007     

Conformational memories with variable bond angles

Conformational Memories (CM) is a simulated annealing/Monte Carlo method that explores peptide and protein dihedral conformational space completely and efficiently, independent of the original conformation. Here we extend the CM method to include the variation of a randomly chosen bond angle, in addition to the standard variation of two or three randomly chosen dihedral angles, in each Monte Carlo trial of the CM exploratory and biased phases. We test the hypothesis that the inclusion of variable bond angles in CM leads to an improved sampling of conformational space. We compare the results with variable bond angles to CM with no bond angle variation for the following systems: (1) the pentapeptide Met-enkephalin, which is a standard test case for conformational search methods; (2) the proline ring pucker in a 17mer model peptide, (Ala)(8)Pro(Ala)(8); and (3) the conformations of the Ser 7.39 chi(1) in transmembrane helix 7 (TMH7) of the cannabinoid CB1 receptor, a 25-residue system. In each case, analysis of the CM results shows that the inclusion of variable bond angles results in sampling of regions of conformational space that are inaccessible to CM calculations with only variable dihedral angles, and/or a shift in conformational populations from those calculated when variable bond angles are not included. The incorporation of variable bond angles leads to an improved sampling of conformational space without loss of efficiency. Our examples show that this improved sampling leads to better exploration of biologically relevant conformations that have been experimentally validated.     

Molecule‐adapted basis sets optimized with a quantum Monte Carlo method

The Monte Carlo simulated annealing method is adapted to optimize correlated Gaussian‐type functions in nonrelativistic molecular environments. Starting from an atom‐centered atomic Gaussian basis set, the uncontracted functions are reoptimized in the molecular environments corresponding to the H₂O, CN^?, N₂, CO, BF, NO⁺, CO₂, and CS systems. These new molecular adapted basis sets are used to calculate total energies, harmonic vibrational frequencies, and equilibrium geometries at a correlated level of theory. The present methodology is a simple and effective way to improve molecular correlated wave functions, without the need to enlarge the molecular basis set. Additionally, this methodology can be used to generate hierarchical sequences of molecular basis sets with increasing size, which are relevant to establish complete basis set limits. © 2014 Wiley Periodicals, Inc.     

Some recent applications of ab initio electronic structure methods to metal,semimetal, and molecular clusters

Recent experimental and theoretical cluster studies are reviewed. Areas of current and developing interest in theoretical and computational chemistry are identified. Some promising methods applied to metal clusters, main group clusters, molecular clusters, spectroscopy, and models of cluster-molecule reactions are indicated. Results of calculations on small hydrogenated lithium clusters and hydrated sodium clusters are discussed in some detail.     

Systematic stepsize variation: Efficient method for searching conformational space of polypeptides

A new and efficient method for overcoming the multiple minima problem of polypeptides, the systematic stepsize variation (SSV) method, is presented. The SSV is based on the assumption that energy barriers can be passed over by sufficiently large rotations about rotatable bonds: randomly chosen dihedral angles are updated starting with a small stepsize (i.e., magnitude of rotation). A new structure is accepted only if it possesses a lower energy than the precedent one. Local minima are passed over by increasing the stepsize systematically. When no new structures are found any longer, the simulation is continued with the starting structure, but other trajectories will be followed due to the random order in updating the torsional angles. First, the method is tested with Met-enkephalin, a peptide with a known global minimum structure; in all runs the latter is found at least once. The global minimum conformations obtained in the simulations show deviations of ±0.0004 kcal/mol from the reference structure and, consequently, are perfectly superposable. For comparison, Metropolis Monte Carlo simulated annealing (MMC-SA) is performed. To estimate the efficiency of the algorithm depending on the complexity of the optimization problem, homopolymers of Ala and Gly of different lengths are simulated, with both the SSV and the MMC-SA method. The comparative simulations clearly reveal the higher efficiency of SSV compared with MMC-SA. © 1998 John Wiley & Sons, Inc. J Comput Chem 19: 1470–1481, 1998     

On the importance of carbohydrate-aromatic interactions for the molecular recognition of oligosaccharides by proteins: NMR studies of the structure and binding affinity of AcAMP2-like peptides with non-natural naphthyl and fluoroaromatic residues

The specific interaction of a variety of modified hevein domains to chitooligosaccharides has been studied by NMR spectroscopy in order to assess the importance of aromatic-carbohydrate interactions for the molecular recognition of neutral sugars. These mutant AcAMP2-like peptides, which have 4-fluoro-phenylalanine, tryptophan, or 2-naphthylalanine at the key interacting positions, have been prepared by solid-phase synthesis. Their three-dimensional structures, when bound to the chitin-derived trisaccharide, have been deduced by NMR spectroscopy. By using DYANA and restrained molecular dynamics simulations with the AMBER 5.0 force field, the three-dimensional structures of the protein-sugar complexes have been obtained. The thermodynamic analysis of the interactions that occur upon complex formation have also been carried out. Regarding binding affinity, the obtained data have permitted the deduction that the larger the aromatic group, the higher the association constant and the binding enthalpy. In all cases, entropy opposes binding. In contrast, deactivation of the aromatic rings by attaching fluorine atoms decreases the binding affinity, with a concomitant decrease in enthalpy. The role of the chemical nature of the aromatic ring for establishing sugar contacts has been thus evaluated.     

Walking freely in the energy and temperature space by the modified replica exchange molecular dynamics method

Replica Exchange Molecular Dynamics (REMD) method is a powerful sampling tool in molecular simulations. Recently, we made a modification to the standard REMD method. It places some inactive replicas at different temperatures as well as the active replicas. The method completely decouples the number of the active replicas and the number of the temperature levels. In this article, we make a further modification to our previous method. It uses the inactive replicas in a different way. The inactive replicas first sample in their own knowledge‐based energy databases and then participate in the replica exchange operations in the REMD simulation. In fact, this method is a hybrid between the standard REMD method and the simulated tempering method. Using different active replicas, one can freely control the calculation quantity and the convergence speed of the simulation. To illustrate the performance of the method, we apply it to some small models. The distribution functions of the replicas in the energy space and temperature space show that the modified REMD method in this work can let the replicas walk freely in both of the two spaces. With the same number of the active replicas, the free energy surface in the simulation converges faster than the standard REMD. © 2016 Wiley Periodicals, Inc.     

Assessing search strategies for flexible docking

We assess the efficiency of molecular dynamics (MD), Monte Carlo (MC), and genetic algorithms (GA) for docking five representative ligand–receptor complexes. All three algorithms employ a modified CHARMM-based energy function. The algorithms are also compared with an established docking algorithm, AutoDock. The receptors are kept rigid while flexibility of ligands is permitted. To test the efficiency of the algorithms, two search spaces are used: an 11-Å-radius sphere and a 2.5-Å-radius sphere, both centered on the active site. We find MD is most efficient in the case of the large search space, and GA outperforms the other methods in the small search space. We also find that MD provides structures that are, on average, lower in energy and closer to the crystallographic conformation. The GA obtains good solutions over the course of the fewest energy evaluations. However, due to the nature of the nonbonded interaction calculations, the GA requires the longest time for a single energy evaluation, which results in a decreased efficiency. The GA and MC search algorithms are implemented in the CHARMM macromolecular package. © 1998 John Wiley & Sons, Inc. J Comput Chem 19: 1623–1631, 1998