Comparative analysis of protein-bound ligand conformations with respect to catalyst's conformational space subsampling algorithms

We examined the quality of Catalyst's conformational model generation algorithm via a large scale study based on the crystal structures of a sample of 510 pharmaceutically relevant protein-ligand complexes extracted from the Protein Data Bank (PDB). Our results show that the tested algorithms implemented within Catalyst are able to produce high quality conformers, which in most of the cases are well suited for in silico drug research. Catalyst-specific settings were analyzed, such as the method used for the conformational model generation (FAST vs BEST) and the maximum number of generated conformers. By setting these options for higher fitting quality, the average RMS values describing the similarity of experimental and simulated conformers were improved from an RMS of 1.06 with max. 50 FAST generated conformers to an RMS of 0.93 with max. 255 BEST generated conformers, which represents an improvement by 12%. Each method provides best fitting conformers with an RMS value<1.50 in more than 80% of all cases. We analyzed the computing time/quality ratio of various conformational model generation settings and examined ligands in high energy conformations. Furthermore, properties of the same ligands in various proteins were investigated, and the fitting qualities of experimental conformations from the PDB and the Cambridge Structural Database (CSD) were compared. One of the most important conclusions of former studies, the fact that bioactive conformers often have energy high above that of global minima, was confirmed.     

High temperature annealed molecular dynamics simulations as a tool for conformational sampling. Application to the bicyclic “222” cryptand

We have performed “High Temperature Annealed Molecular Dynamics Simulations” (HTAMDS) on the bicyclic 222 cryptand, and on model M⁺/222 cryptates with different representations of M⁺. The analysis of four sets of 500 structures allows assessment of the ability of HTAMDS to: (1) interconvert experimentally known conformers starting from one of them, (2) locate the energy minima, (3) generate new conformers of low energy, and (4) account for the average structure observed on the NMR time scale. In view of the ionophoric behavior of 222, structures are also analyzed in terms of the “in/out” orientation of the binding sites. It is found that simulations on the free molecule, although widely sampling the conformational space, do not give structures adequate for cation inclusion. They generate however the lowest energy structure known experimentally and other new closely related ones. Inclusion of the substrate in the simulation (either as a purely electrostatic “driver,” or as a charged sphere) is required to generate conformations found in several complexes. These results suggest that in the field of drug design, conformation of drugs suitable for binding to a given receptor may not be found when the simulations are performed on the isolated drug or substrate.     

Searching and optimizing structure ensembles for complex flexible sugars

NMR restrictions are suitable to specify the geometry of a molecule when a single well-defined global free energy minimum exists that is significantly lower than other local minima. Carbohydrates are quite flexible, and therefore, NMR observables do not always correlate with a single conformer but instead with an ensemble of low free energy conformers that can be accessed by thermal fluctuations. In this communication, we describe a novel procedure to identify and weight the contribution to the ensemble of local minima conformers based on comparison to residual dipolar couplings (RDCs) or other NMR observables, such as scalar couplings. A genetic algorithm is implemented to globally minimize the R factor comparing calculated RDCs to experiment. This is done by optimizing the weights of different conformers derived from the exhaustive local minima conformational search program, fast sugar structure prediction software (FSPS). We apply this framework to six human milk sugars, LND-1, LNF-1, LNF-2, LNF-3, LNnT, and LNT, and are able to determine corresponding population weights for the ensemble of conformers. Interestingly, our results indicate that in all cases the RDCs can be well represented by only a few most important conformers. This confirms that several, but not all of the glycosidic linkages in histo-blood group "epitopes" are quite rigid.     

A density functional study of the structure and potential energy surface of dibenzenechromium and its alkylated derivatives

Density functional calculations were performed with the nonempirically constructed PBE functional and TZ2p basis set to study the potential energy surfaces of several dibenzenechromium derivatives containing methy1 and tert-buty1 groups in aromatic ligands. The method was shown to correctly describe the structure and intramolecular dynamics of bis-arenechromium complexes. It allows barriers to rotation about the metal-ligand bond to be calculated and detailed information to be obtained about the shape of the potential curves that correspond to this rotation. The potential energy surfaces of polymethylated dibenzenechromium derivatives contained several minima close in energy; these minima largely corresponded to eclipsed conformers. The potential curves for rotation about the metal-ligand bond in dibenzenechromium and its methylated derivatives reveal the presence of staggered conformers as transition states. Irrespective of the number of methyl groups present in both ligands, we observed no substantial changes in barrier heights, which were about 1 kcal/mol. Conversely, the introduction of bulky tert-butyl groups destabilized eclipsed conformations. Ensuing steric strain caused substantial out-of-ring-plane displacements of the C_ipso-t-Bu bonds away from the chromium atom, distorted aromatic ring planarity, and substantially increased barriers to rotation, to 8–10 kcal/mol.     

Ab initio calculations on low-energy conformers of alpha-cyclodextrin

In this article we carried out a comprehensive investigation of true minima on the potential energy surface (PES) for the alpha-cyclodextrin molecule using ab initio Hartree-Fock (HF) and density functional theory (DFT) quantum chemical methods, employing basis sets ranging from 6-31 G(d,p) to 6-311++G(2d,2p) triple-zeta quality. Thermodynamic quantities and the solvent effect were evaluated at the DFT level of theory. We believe that the most relevant conformers present on the multidimensional PES were sampled in our work, using an adequate treatment of electron correlation effects to describe the intramolecular hydrogen bonds that are present in cyclodextrin species. We present new structures not reported so far and discuss, in detail, the relevance of the DFT gas-phase equilibrium structures for the experimental and theoretical studies involving cyclodextrins and corresponding inclusion complexes, in the condensed phase. In addition, among the various true minimum energy structures located on the DFT PES, the preferred structures in the gas phase and aqueous media, needed to be used as representative minima on the PES in further studies involving the interaction of alpha-cyclodextrin with other species, were unambiguously identified.     

An investigation of the adiabatic potential energy surface of a Pseudo-Jahn-Teller complex using X-ray diffraction

X-ray diffraction combined with theoretical results obtained for isolated molecules can reveal many topologic features of the adiabatic potential energy surface (APS) of pseudo-Jahn-Teller Cu(II) complexes. Single-crystal x-ray results of acetatobis (1,10-phenanthroline)copper(II) perchlorate from 10 to 351 K are presented. The observed Cu-O bond length changes of 0.26 and 0.38 Å are rationalized in terms of a one-dimensional APS with two nonequivalent minima. The Silver and Getz model of vibronic coupling applied to the temperature variation of the Cu-O bond lengths is not obeyed. The actual structures of the conformers associated with the energy minima were determined from the x-ray data. A means of estimating the pseudo-Jahn-Teller radius using mean-square vibrational amplitudes along the Cu-O bonds is also discussed.     

Ab initio studies of conformational properties of dimethyl diphosphate dianion and its complex with magnesium

The conformational properties of the diphosphate linkage have been studied with ab initio methods using the dimethyl diphosphate dianion (1) and magnesium dimethyl diphosphate (2) as models. The ab initio energy and geometry of the conformers around the P-O bonds have been determined at the self-consistent-field (SCF) using the 6-31G* and the tzp basis sets; whereas, the 6-31G* basis set alone has been used for 2. In addition, the adiabatic connection method (ACM) of density functional theory (DFT) using the dzvp basis set has been employed for 1. The optimization of all possible staggered conformers assumed for the four P-O bonds, led to nine minima for 1. In agreement with the general anomeric effect, the sc conformation about the P-O bonds is clearly preferred over the ap one. Vibrational frequencies were calculated at the SCF level using the 6-31G* basis set and used to evaluate zero-point energies, thermal energies, and entropies for all minima of 1. The effect of zero-point energies and thermal energies is quite small. However, the effect of entropies, mainly resulting from a multiplicity contribution, changes the stability of the conformers. For each minimum of 1, up to six different arrangements of the Mg²⁺ were used to determine minima of 2. This procedure led to 21 distinct minima. The presence of the magnesium counter-ion appeared to completely change the structure and relative energy of the conformers. The preferred structures of the complex exhibit the (sc, ap) orientation around the two central P–O bonds and an arrangement in which the magnesium cation is coordinated by three phosphoryl oxygen atoms. The results of this work clearly demonstrate that interactions with the metal counter-ion can induce conformational changes in the overall 3D-shape adopted by molecules containing diphosphate linkages. The PM3 and MNDO quantum semi-empirical methods and molecular mechanics methods using the CVFF force field were tested and large differences in the minimum structures, as well as in the conformational energies between these and ab initio methods, are discussed.     

Conformational analysis of macrocycles: comparing general and specialized methods

Macrocycles represent an important class of medicinally relevant small molecules due to their interesting biological properties. Therefore, a firm understanding of their conformational preferences is important for drug design. Given the importance of macrocycle-protein modelling in drug discovery, we envisaged that a systematic study of both classical and recent specialized methods would provide guidance for other practitioners within the field. In this study we compare the performance of the general, well established conformational analysis methods Monte Carlo Multiple Minimum (MCMM) and Mixed Torsional/Low-Mode sampling (MTLMOD) with two more recent and specialized macrocycle sampling techniques: MacroModel macrocycle Baseline Search (MD/LLMOD) and Prime macrocycle conformational sampling (PRIME-MCS). Using macrocycles extracted from 44 macrocycle-protein X-ray crystallography complexes, we evaluated each method based on their ability to (i) generate unique conformers, (ii) generate unique macrocycle ring conformations, (iii) identify the global energy minimum, (iv) identify conformers similar to the X-ray ligand conformation after Protein Preparation Wizard treatment (X-rayppw), and (v) to the X-rayppw ring conformation. Computational speed was also considered. In addition, conformational coverage, as defined by the number of conformations identified, was studied. In order to study the relative energies of the bioactive conformations, the energy differences between the global energy minima and the energy minimized X-rayppw structures and, the global energy minima and the MCMM-Exhaustive (1,000,000 search steps) generated conformers closest to the X-rayppw structure, were calculated and analysed. All searches were performed using relatively short run times (10,000 steps for MCMM, MTLMOD and MD/LLMOD). To assess the performance of the methods, they were compared to an exhaustive MCMM search using 1,000,000 search steps for each of the 44 macrocycles (requiring ca 200 times more CPU time). Prior to our analysis, we also investigated if the general search methods MCMM and MTLMOD could also be optimized for macrocycle conformational sampling. Taken together, our work concludes that the more general methods can be optimized for macrocycle modelling by slightly adjusting the settings around the ring closure bond. In most cases, MCMM and MTLMOD with either standard or enhanced settings performed well in comparison to the more specialized macrocycle sampling methods MD/LLMOD and PRIME-MCS. When using enhanced settings for MCMM and MTLMOD, the X-rayppw conformation was regenerated with the greatest accuracy. The, MD/LLMOD emerged as the most efficient method for generating the global energy minima.     

Potential-energy and free-energy surfaces of glycyl-phenylalanyl-alanine (GFA) tripeptide: experiment and theory

The free-energy surface (FES) of glycyl-phenylalanyl-alanine (GFA) tripeptide was explored by molecular dynamics (MD) simulations in combination with high-level correlated ab initio quantum chemical calculations and metadynamics. Both the MD and metadynamics employed the tight-binding DFT-D method instead of the AMBER force field, which yielded inaccurate results. We classified the minima localised in the FESs as follows: a) the backbone-conformational arrangement; and b) the existence of a COOH...OC intramolecular H-bond (families CO(2)H(free) and CO(2)H(bonded)). Comparison with experimental results showed that the most stable minima in the FES correspond to the experimentally observed structures. Remarkably, however, we did not observe experimentally the CO(2)H(bonded) family (also predicted by metadynamics), although its stability is comparable to that of the CO(2)H(free) structures. This fact was explained by the former's short excited-state lifetime. We also carried out ab initio calculations using DFT-D and the M06-2X functional. The importance of the dispersion energy in stabilising peptide conformers is well reflected by our pioneer analysis using the DFT-SAPT method to explore the nature of the backbone/side-chain interactions.     

Conformation resolved induced infrared activity: trans- and cis-formic acid isolated in solid molecular hydrogen

We report combined experimental and theoretical studies of infrared absorptions induced in solid molecular hydrogen by different conformers of formic acid (HCOOH, FA). FTIR spectra recorded in the H(2) fundamental region (4120-4160 cm(-1)) reveal a number of relatively strong trans-FA induced Q-branch absorptions that are assigned by studying both FA-doped parahydrogen (pH(2)) and normal hydrogen (nH(2)) samples. The induced H(2) absorptions are also studied for HCOOD doped nH(2) crystals for both the trans and cis conformers that show resolvable differences. Samples containing >90% of the higher energy cis-HCOOD conformer are produced by in situ IR pumping of the OD stretching overtone of trans-HCOOD using narrow-band IR light. Minimum energy structures for 1:1 complexes of H(2) and FA are determined using ab initio methods. The measured differences in the cis- versus trans-HCOOD induced spectra are in qualitative agreement with the frequencies and intensities calculated for the identified cluster structures as discussed in terms of the model of specific interactions.     

Searching for conformers of nine- to twelve-ring hydrocarbons on the MM2 and MM3 energy surfaces: Stochastic search for interconversion pathways

The stochastic search method was employed to find as many conformers on the MM2 and MM3 energy surfaces as possible for cyclic saturated hydrocarbons with ring sizes from 9 through 12. The number found was 8 MM2 (8 MM3) for 9 rings, 18 MM2 (16 MM3) for 10 rings, 40 MM2 (29 MM3) for 11 rings, and 111 MM2 (90 MM3) for 12 rings. A measure of similarity between pairs of conformers of a compound, called conformational distance, is described. It was used to correlate similar MM2 and MM3 conformers. It was discovered that some conformers on each energy surface are not close to minima on the other surface in rings larger than 9. On refinement with the other optimizer, they changed considerably—going downhill to other previously found minima on the other energy surface or (in a few cases) going to minima which had not been found by direct searches. Conformational distance was also employed as an indication of which pairs of MM2 (or MM3) conformers are likely to interconvert rapidly. A new stochastic procedure of using small kicks was used to search for the most likely interconversion processes among the conformers. There is fairly good agreement between the most facile pathways located by it and unusually short conformational distances. Several additional 12-ring conformers (not found with previous methods) were located through application of this small kick procedure.     

Monte Carlo temperature basin paving with effective fragment potential: an efficient and fast method for finding low-energy structures of water clusters (H2O)20 and (H2O)25

Determining low-energy structures of large water clusters is a challenge for any optimization algorithm. In this work, we have developed a new Monte Carlo (MC)-based method, temperature basin paving (TBP), which is related to the well-known basin hopping method. In the TBP method, the Boltzmann weight factor used in MC methods is dynamically modified based on the history of the simulation. The states that are visited more are given a lower probability by increasing their temperatures and vice versa. This allows faster escapes from the states frequently visited in the simulation. We have used the TBP method to find a large number of low-energy minima of water clusters of size 20 and 25. We have found structures energetically same to the global minimum structures known for these two clusters. We have compared the efficiency of this method to the basin-hopping method and found that it can locate the minima faster. Statistical efficiency of the new method has been investigated by running a large number of trajectories. The new method can locate low-energy structures of both the clusters faster than some of the reported algorithms for water clusters and can switch between high energy and low-energy structures multiple times in a simulation illustrating its efficiency. The large number of minima obtained from the simulations is used to get both general and specific features of the minima. The distribution of minima for these two clusters based on the similarity of their oxygen frames shows that the (H(2)O)(20) can have different variety of structures, but for (H(2)O)(25), low-energy structures are mostly cagelike. Several (H(2)O)(25) structures are found with similar energy but with different cage architectures. Noncage structures of (H(2)O)(25) are also found but they are 6-7 kcal/mol higher in energy from the global minimum. The TBP method is likely to play an important role for exploring the complex energy landscape of large molecules.     

Ab initio study of ascorbic acid conformations

The STO-3G optimized structures of nine different staggered conformers of ascorbic acid are presented. The largest energy difference between the nine local minima is 5.1 kcal/mol. Comparison of the relative energies of the fully optimized structures of ascorbic acid conformers with those of nonoptimized conformers shows that full optimization is essential to obtain meaningful results. However, optimization of the ring structure is almost independent of optimization of the side-chain structure. One of the STO-3G optimized gas phase conformers is very close to the X-ray structure of the crystal.     

Potential energy and free energy landscapes

Familiar concepts for small molecules may require reinterpretation for larger systems. For example, rearrangements between geometrical isomers are usually considered in terms of transitions between the corresponding local minima on the underlying potential energy surface, V. However, transitions between bulk phases such as solid and liquid, or between the denatured and native states of a protein, are normally addressed in terms of free energy minima. To reestablish a connection with the potential energy surface we must think in terms of representative samples of local minima of V, from which a free energy surface is projected by averaging over most of the coordinates. The present contribution outlines how this connection can be developed into a tool for quantitative calculations. In particular, stepping between the local minima of V provides powerful methods for locating the global potential energy minimum, and for calculating global thermodynamic properties. When the transition states that link local minima are also sampled we can exploit statistical rate theory to obtain insight into global dynamics and rare events. Visualizing the potential energy landscape helps to explain how the network of local minima and transition states determines properties such as heat capacity features, which signify transitions between free energy minima. The organization of the landscape also reveals how certain systems can reliably locate particular structures on the experimental time scale from among an exponentially large number of local minima. Such directed searches not only enable proteins to overcome Levinthal's paradox but may also underlie the formation of "magic numbers" in molecular beams, the self-assembly of macromolecular structures, and crystallization.     

Improved mapping of protein binding sites

Computational mapping methods place molecular probes – small molecules or functional groups – on a protein surface in order to identify the most favorable binding positions by calculating an interaction potential. Mapping is an important step in a number of flexible docking and drug design algorithms. We have developed improved algorithms for mapping protein surfaces using small organic molecules as molecular probes. The calculations reproduce the binding of eight organic solvents to lysozyme as observed by NMR, as well as the binding of four solvents to thermolysin, in good agreement with x-ray data. Application to protein tyrosine phosphatase 1B shows that the information provided by the mapping can be very useful for drug design. We also studied why the organic solvents bind in the active site of proteins, in spite of the availability of alternative pockets that can very tightly accommodate some of the probes. A possible explanation is that the binding in the relatively large active site retains a number of rotational states, and hence leads to smaller entropy loss than the binding elsewhere else. Indeed, the mapping reveals that the clusters of the ligand molecules in the protein's active site contain different rotational-translational conformers, which represent different local minima of the free energy surface. In order to study the transitions between different conformers, reaction path and molecular dynamics calculations were performed. Results show that most of the rotational states are separated by low free energy barriers at the experimental temperature, and hence the entropy of binding in the active site is expected to be high.     

Theoretical study of hydrated copper(II) interactions with guanine: a computational density functional theory study

Optimization of the hydrated Cu(II)(N7-guanine) structures revealed a number of minima on the potential energy surface. For selected structures, energy decompositions together with the determination of electronic properties (partial charges and electron spin densities) were performed. In the complexes of guanine with the bare copper cation and that with the monoaqua ligated cation, an electron transfer from guanine to Cu(II) was observed, resulting in a Cu(I)-guanine(+) type of complex. Conformers with two aqua ligands are borderline systems characterized by a Cu partial charge of +0.7e and a similar value of the spin density (0.6e) localized on guanine. When tetracoordination of copper was achieved, only then the prevailing electron spin density is unambiguously localized on copper. The energetic preference of diaqua-Cu-(N7,O6-guanine) over triaqua-Cu-(N7-guanine) was found for the four-coordinate structures. However, the energy difference between these two conformations decreases with the number of water molecules present in the systems, and in complexes with five water molecules this preference is preserved only at DeltaG level where thermal and entropy terms are included.     

Substituent effects on the photophysical properties of 2,9-substituted phenanthroline copper(I) complexes: a theoretical investigation

The electronic and nuclear structures of a series of [Cu(2,9-(X)₂-phen)₂]⁺ copper(I) complexes (phen=1,10-phenanthroline; X=H, F, Cl, Br, I, Me, CN) in their ground and excited states are investigated by means of density functional theory (DFT) and time-dependent (TD-DFT) methods. Subsequent Born-Oppenheimer molecular dynamics is used for exploring the T₁ potential energy surface (PES). The T₁ and S₁ energy profiles, which connect the degenerate minima induced by ligand flattening and Cu−N bond symmetry breaking when exciting the molecule are calculated as well as transition state (TS) structures and related energy barriers. Three nuclear motions drive the photophysics, namely the coordination sphere asymmetric breathing, the well-documented pseudo Jahn-Teller (PJT) distortion and the bending of the phen ligands. This theoretical study reveals the limit of the static picture based on potential energy surfaces minima and transition states for interpreting the luminescent and TADF properties of this class of molecules. Whereas minor asymmetric Cu−N bonds breathing accompanies the metal-to-ligand-charge-transfer re-localization over one or the other phen ligand, the three nuclear movements participate to the flattening of the electronically excited complexes. This leads to negligible energy barriers whatever the ligand X for the first process and significant ligand dependent energy barriers for the formation of the flattened conformers. Born-Oppenheimer (BO) dynamics simulation of the structural evolution on the T₁ PES over 11 ps at 300 K confirms the fast backwards and forwards motion of the phenanthroline within 200–300 fs period and corroborates the presence of metastable C₂ structures.     

Gas Phase Conformations of Tetrapeptide Glycine-Phenylalanine-Glycine-Glycine

Systematic search of the potential energy surface of tetrapeptide glycine-phenylalanine-glycine-glycine (GFGG) in gas phase is conducted by a combination of PM3, HF and BHandHLYP methods. The conformational search method is described in detail. The rela-tive electronic energies, zero point vibrational energies, dipole moments, rotational constants, vertical ionization energies and the temperature dependent conformational distributions for a number of important conformers are obtained. The structural characteristics of these conformers are analyzed and it is found that the entropic effect is a dominating factor in determining the relative stabilities of the conformers. The measurements of dipole moments and some characteristic IR mode are shown to be effective approaches to verify the theoreti-cal prediction. The structures of the low energy GFGG conformers are also analyzed in their connection with the secondary structures of proteins. Similarity between the local structures of low energy GFGG conformers and the α-helix is discussed and many β- and γ-turn local structures in GFGG conformers are found.