Conformational searching methods for small molecules. II. Genetic algorithm approach

We demonstrate the use of a genetic algorithm (GA) search procedure for finding low-energy conformations of small to medium organic molecules (1–12 rotatable bonds). GAS are in a class of biologically motivated optimization methods that evolve a population of individuals where individuals who are more “fit” have a higher probability of surviving into subsequent generations. Here, an individual is a conformation of a given molecule and the fitness is the molecule's conformational energy. In the course of a simulated evolution, the population produces conformations having increasingly lower energy. We test the GA method on a suite of 72 molecules and compare the performance against the CSEARCH algorithm in Sybyl. For molecules with more than eight rotatable bonds, the GA method is more efficient computationally and as the number of rotatable bonds increases the relative efficiency of the GA method grows. The GA method also found energies equal to or lower than the energy of the relaxed crystal structure in the large majority of cases. © John Wiley & Sons, Inc.     

The development of a simple empirical scoring function to estimate the binding constant for a protein-ligand complex of known three-dimensional structure

Summary A new simple empirical function has been developed that estimates the free energy of binding for a given protein-ligand complex of known 3D structure. The function takes into account hydrogen bonds, ionic interactions, the lipophilic protein-ligand contact surface and the number of rotatable bonds in the ligand. The dataset for the calibration of the function consists of 45 protein-ligand complexes. The new energy function reproduces the binding constants (ranging from 2.5·10^-2 to 4·10^-14 M, corresponding to binding energies between -9 and -76 kJ/mol) of the dataset with a standard deviation of 7.9 kJ/mol, corresponding to 1.4 orders of magnitude in binding affinity. The individual contributions to protein-ligand binding obtained from the scoring function are: ideal neutral hydrogen bond: -4.7 kJ/mol; ideal ionic interaction: -8.3 kJ/mol; lipophilic contact: -0.17 kJ/mol Ų; one rotatable bond in the ligand: +1.4 kJ/mol. The function also contains a constant contribution (+5.4 kJ/mol) which may be rationalized as loss of translational and rotational entropy. The function can be evaluated very fast and is therefore also suitable for application in a 3D database search or de novo ligand design program such as LUDI.     

Generating conformer ensembles using a multiobjective genetic algorithm

The task of generating a nonredundant set of low-energy conformations for small molecules is of fundamental importance for many molecular modeling and drug-design methodologies. Several approaches to conformer generation have been published. Exhaustive searches suffer from the exponential growth of the search space with increasing degrees of conformational freedom (number of rotatable bonds). Stochastic algorithms do not suffer as much from the exponential increase of search space and provide a good coverage of the energy minima. Here, the use of a multiobjective genetic algorithm in the generation of conformer ensembles is investigated. Distance geometry is used to generate an initial conformer, which is then subject to geometric modifications encoded by the individuals of the genetic algorithm. The geometric modifications apply to torsion angles about rotatable bonds, stereochemistry of double bonds and tetrahedral chiral centers, and ring conformations. The geometric diversity of the evolving conformer ensemble is preserved by a fitness-sharing mechanism based on the root-mean-square distance of the atomic coordinates. Molecular symmetry is taken into account in the distance calculation. The geometric modifications introduce strain into the structures. The strain is relaxed using an MMFF94-like force field in a postprocessing step that also removes conformational duplicates and structures whose strain energy remains above a predefined window from the minimum energy value found in the set. The implementation, called Balloon, is available free of charge on the Internet ( http://www.abo.fi/~mivainio/balloon/).     

Molecular properties from conformational ensembles. 1. Dipole moments of molecules with multiple internal rotations

We present a novel method for constructing the stable conformational space of small molecules with many rotatable bonds that uses our iterative stochastic elimination (ISE) algorithm, a robust stochastic search method capable of finding ensembles of best solutions for large combinatorial problems. To validate the method, we show that ISE reproduces the best conformers found in a fully exhaustive search, as well as compare computed dipole moments to experimental values, based on molecular ensembles and their Boltzmann distributions. Results were also compared to the alternative molecular dynamics and simulated annealing methods. Our results clarify that many low energy conformations may be required to reproduce molecular properties, while single low energy conformers or ensembles of low energy conformers cannot account for the experimental properties of flexible molecules. Whereas ISE well reproduces conformations that are not separated by very large energy barriers, it has not been successful in reproducing conformations of strained molecules.     

Conformational energy downward driver (CEDD): Characterization and calibration of the method

Summary A method has been developed that allows one to drive a molecule to conformations of lowest energy given the starting conformation, the identity of the rotatable bonds and the step size. This method has proved useful in our hands in the drug design arena where it is frequently more important to get low-energy conformers of a molecule that match some other (e.g. pharmacophoric or enzyme pocket) requirements than to exhaustively enumerate all possible low-energy conformations for each of the molecules to be studied. The method has been shown to work in the test cases studied to date. Furthermore, so far it has been shown to be sufficiently fast to be used for molecules containing up to 70 rotatable bonds.     

Understanding Ramachandran plot for dipeptide: A density functional theory and information-theoretic approach study

The Ramachandran plot in biochemistry visualizes the propensity of amino acid residues to assume different secondary structures. It can be obtained by different approaches such as bioinformatics and molecular dynamics. In this work, we systematically scan the two single bonds around the α-carbon of a peptide to generate the three-dimensional potential energy surface and then project it onto the coordinate plane to generate the Ramachandran plot. Using the two schemes of the total energy partition and information-theoretic approach in density functional theory, we analyze the plot and find out that the dominant contributor to determine the local minima of secondary structures, especially β-sheets, in the Ramachandran plot is the electrostatic interaction, whereas steric and exchange-correlation contributions play minor yet indispensable roles, especially for determining α-helixes. As the generalization of our prior studies for systems with one rotatable bond, our current results confirm what we yielded before from the viewpoint of energetic contributions, and at meanwhile provide an in-depth understanding about the nature and origin of secondary structure propensities from the Ramachandran plot.     

Transition metal-carbon complexes. A theoretical study

The equilibrium geometries and bond dissociation energies of 16VE and 18VE complexes of ruthenium and iron with a naked carbon ligand are reported using density functional theory at the BP86/TZ2P level. Bond energies were also calculated at CCSD(T) using TZ2P quality basis sets. The calculations of [Cl2(PMe3)2Ru(C)] (1Ru), [Cl2(PMe3)2Fe(C)] (1Fe), [(CO)2(PMe3)2Ru(C)] (2Ru), [(CO)2(PMe3)2Fe(C)] (2Fe), [(CO)4Ru(C)] (3Ru), and [(CO)4Fe(C)] (3Fe) show that 1Ru has a very strong Ru-C bond which is stronger than the Fe-C bond in 1Fe. The metal-carbon bonds in the 18VE complexes 2Ru-3Fe are weaker than those in the 16VE species. Calculations of the related carbonyl complexes [(PMe3)2Cl2Ru(CO)] (4Ru), [(PMe3)2Cl2Fe(CO)] (4Fe), [(PMe3)2Ru(CO)3] (5Ru), [(PMe3)2Fe(CO)3] (5Fe), [Ru(CO)5] (6Ru), and [Fe(CO)5] (6Fe) show that the metal-CO bonds are much weaker than the metal-C bonds. The 18VE iron complexes have a larger BDE than the 18VE ruthenium complexes, while the opposite trend is calculated for the 16VE compounds. Charge and energy decomposition analyses (EDA) have been carried out for the calculated compounds. The Ru-C and Fe-C bonds in 1Ru and 1Fe are best described in terms of two electron-sharing bonds with sigma and pi symmetry and one donor-acceptor pi bond. The bonding situation in the 18 VE complexes 2Ru-3Fe is better described in terms of closed shell donor-acceptor interactions in accordance with the Dewar-Chatt-Duncanson model. The bonding analysis clearly shows that the 16VE carbon complexes 1Ru and 1Fe are much more strongly stabilized by metal-C sigma interactions than the 18VE complexes which is probably the reason why the substituted homologue of 1Ru could become isolated. The EDA calculations show that the nature of the TM-C and TM-CO binding interactions resembles each other. The absolute values for the energy terms which contribute to Delta(Eint) are much larger for the carbon complexes than for the carbonyl complexes, but the relative strengths of the energy terms are not very different from each other. The pi bonding contribution to the orbital interactions in the carbon complexes is always stronger than sigma bonding. There is no particular bonding component which is responsible for the reversal of the relative bond dissociation energies of the Ru and Fe complexes when one goes from the 16VE complexes to the 18VE species. That the 18 VE compounds have longer and weaker TM-C and TM-CO bonds than the respective 16 VE compounds holds for all complexes. This is because the LUMO in the 16 VE species is a sigma-antibonding orbital which becomes occupied in the 18 VE species.     

Stabilities of Hydrogen-Bonded Supramolecular Complexes with Various Numbers of Single Bonds: Attempts To Quantify a Dogma in Host–Guest Chemistry

The disadvantage of flexible bonds in supramolecular host–guest complexes is much smaller than usually assumed. In the association of dicarboxlic acids with diamides (shown on the right), freely rotatable single bonds have only a minor disadvantageous influence on the free energy of complexation ΔG_cplx. From the linear correlation between ΔG_cplx and the number of single bonds n, a decrease in the free energy of complexation of only 1.3 kJ mol⁻¹ per single bond was determined.     

The assembly-inducing laulimalide/peloruside a binding site on tubulin: molecular modeling and biochemical studies with [³H]peloruside A

We used synthetic peloruside A for the commercial preparation of [3H]peloruside A. The radiolabeled compound bound to preformed tubulin polymer in amounts stoichiometric with the polymer's tubulin content, with an apparent K(d) value of 0.35 μM. A less active peloruside A analogue, (11-R)-peloruside A and laulimalide acted as competitive inhibitors of the binding of the [3H]peloruside A, with apparent K(i) values of 9.3 and 0.25 μM, respectively. Paclitaxel, epothilone B, and discodermolide had essentially no ability to inhibit [3H]peloruside A binding, confirming that these compounds bind to a different site on tubulin polymer. We modeled both laulimalide and peloruside A into the binding site on β-tubulin that was identified by Huzil et al. (J. Mol. Biol. 2008, 378, 1016-1030), but our model provides a more reasonable structural basis for the protein-ligand interaction. There is a more complete desolvation of the peloruside A ligand and a greater array of favorable hydrophobic and electrostatic interactions exhibited by peloruside A at its β-tubulin binding site. In addition, the protein architecture in our peloruside A binding model was suitable for binding laulimalide. With the generation of both laulimalide and peloruside A binding models, it was possible to delineate the structural basis for the greater activity of laulimalide relative to peloruside A and to rationalize the known structure-activity relationship data for both compounds.     

EPR Study of Electronic Structure of [CoF<Subscript>6</Subscript>]<Superscript>3−</Superscript>and B<Subscript>18</Subscript>N<Subscript>18</Subscript> Nano Ring Field Effects on Octahedral Complex

Density functional theory calculations (DFT), as well as hybrid methods (B3LYP) for B₁₈N₁₈-[CoF₆]³⁻ complex have been carried out to study the non-bonded interaction. The geometry of the B₁₈N₁₈ has been optimized at B3LYP method with EPR-II basis set and geometry of the [CoF₆]³⁻ have been optimized at B3LYP method with Def2-TZVP basis set and Stuttgart RSC 1997 Effective Core Potential. The electromagnetic interactions of the [CoF₆]³⁻ molecule embedded in the B₁₈N₁₈ Nano ring have been investigated at B3LYP and total atomic charges, spin densities, dipole moment and isotropic Fermi coupling constants parameters in different loops and bonds of the B₁₈N₁₈-[CoF₆]³⁻ system have been calculated. Also NBO analysis such as electronic delocalization between donor and acceptor bonds has been studied by DFT method. Then we have been investigated the lowest unoccupied molecular orbital (LUMO) and the highest occupied molecular orbital (HOMO) for the lowest energy have been derived to estimate the structural stability of the B₁₈N₁₈-[CoF₆]³⁻ system, and the coefficients of s, p and d orbitals of Co-F bonds involved in B₁₈N₁₈-[CoF₆]³⁻.Thus, hybridization of Co and F atoms can be distinguished based on these NBO data. The Gaussian quantum chemistry package is used for all calculations.     

Hydrogen bonds and ionic interactions in Guanidine/Guanidinium complexes: a computational case study

It is frequently said that hydrogen bonds (HBs) are enhanced by ionic interactions and in this article we intend to determine the degree at which this reinforcement happens. Considering our interest in the Guanidine(neutral)/Guanidinium(cation) system and its particular nature, all the possible 1:1 complexes with the Chloride(anion)/Hydrochloric acid(neutral) system have been studied at different levels of computation (B3LYP with 6-31+G* and TZVP basis sets; MP2 with 6-31+G*, 6-311++G** and aug-cc-pVDZ basis sets; CBS-QB3 and G3MP2). The nature of these interactions established in all the systems and, when possible, at all the levels of computation used in this study, has been analyzed using Atoms in Molecules and Natural Bond Orbital methodologies. By examining the interaction energy, the electron density at the bond critical bonds, the atomic energy, the charge transfer, the orbital energy, and the deformation energy we can conclude that HBs are stronger when the ionic interaction is stronger. Thus, both interactions do not work in an independent manner but one reinforces the other to different degrees depending on the nature of the charges present. Several correlations with the interaction energy have been found and a partition of the contributions of both the HB and ionic forces to the total interactions is proposed.     

Dynamics of hyperthermal collisions of O(3P) with CO

The dynamics of O(3P) + CO collisions at a hyperthermal collision energy near 80 kcal mol-1 have been studied with a crossed molecular beams experiment and with quasi-classical trajectory calculations on computed potential energy surfaces. In the experiment, a rotatable mass spectrometer detector was used to monitor inelastically and reactively scattered products as a function of velocity and scattering angle. From these data, center-of-mass (c.m.) translational energy and angular distributions were derived for the inelastic and reactive channels. Isotopically labeled C18O was used to distinguish the reactive channel (16O + C18O 16OC + 18O) from the inelastic channel (16O + C18O 16O + C18O). The reactive 16OC molecules scattered predominantly in the forward direction, i.e., in the same direction as the velocity vector of the reagent O atoms in the c.m. frame. The c.m. translational energy distribution of the reactively scattered 16OC and 18O was very broad, indicating that 16OC is formed with a wide range of internal energies, with an average internal excitation of approximately 40% of the available energy. The c.m. translational energy distribution of the inelastically scattered C18O and 16O products indicated that an average of 15% of the collision energy went into internal excitation of C18O, although a small fraction of the collisions transferred nearly all the collision energy into internal excitation of C18O. The theoretical calculations, which extend previously published results on this system, predict c.m. translational energy and angular distributions that are in near quantitative agreement with the experimentally derived distributions. The theoretical calculations, thus validated by the experimental results, have been used to derive internal state distributions of scattered CO products and to probe in detail the interactions that lead to the observed dynamical behavior.     

Pharmacophore mapping in the laulimalide series: total synthesis of a vinylogue for a late-stage metathesis diversification strategy

An efficient synthesis of the macrocyclic core of laulimalide with a pendant vinyl group at C20 is described, allowing for late-stage introduction of various side chains through a selective and efficient cross metathesis diversification step. Representative analogues reported herein are the first to contain modifications to only the side chain dihydropyran of laulimalide and des-epoxy laulimalide. This step-economical strategy enables the rapid synthesis of new analogues using alkenes as an inexpensive, abundantly available diversification feedstock.     

Treating entropy and conformational changes in implicit solvent simulations of small molecules

Implicit solvent models are increasingly popular for estimating aqueous solvation (hydration) free energies in molecular simulations and other applications. In many cases, parameters for these models are derived to reproduce experimental values for small molecule hydration free energies. Often, these hydration free energies are computed for a single solute conformation, neglecting solute conformational changes upon solvation. Here, we incorporate these effects using alchemical free energy methods. We find significant errors when hydration free energies are estimated using only a single solute conformation, even for relatively small, simple, rigid solutes. For example, we find conformational entropy (TDeltaS) changes of up to 2.3 kcal/mol upon hydration. Interestingly, these changes in conformational entropy correlate poorly (R2 = 0.03) with the number of rotatable bonds. The present study illustrates that implicit solvent modeling can be improved by eliminating the approximation that solutes are rigid.     

Self-assembly of doublets from flattened polymer colloids

Bottom-up fabrication methods are used to assemble strong yet flexible colloidal doublets. Part of a spherical particle is flattened, increasing the effective interaction area with another particle having a flat region. In the presence of a moderate ionic strength, the flat region on one particle will preferentially "bond" to a flat region on another particle in a deep (≥10 kT) secondary energy minimum. No external field is applied during the assembly process. Under the right conditions, the flat-flat bonding strength is ≥10× that of a sphere-sphere interaction. Not only can flat-flat bonds be quite strong, but they are expected to remain freely rotatable and flexible, with negligible energy barriers for rotation because particles reside in a deep secondary energy minimum with a ~20-30 nm layer of fluid between the ~1 μm radius particles. We present a controlled technique to flatten the particles at room temperature, the modeling of the interparticle forces for flattened spheres, and the experimental data for the self-assembly of flat-flat doublets.     

Hydronium ion complex of 18-crown-6: where are the protons? A density functional study of static and dynamic properties

A quantum-chemical study employing the BLYP density functional is reported for the complex of H3O+ with 18-crown-6. According to a Car-Parrinello molecular dynamics (CPMD) study at 340 K, the complex is quite flexible, and is characterized by three quasi-linear (two-center) hydrogen-bond interactions for most of the time. On a time scale of 10 ps, frequent inversions of H3O+ are observed, as well as two 120 degrees rotations switching the hydrogen bonds from one set of crown-ether O atoms to the other. These results are consistent with density-functional studies of stationary points on the potential energy surface, which show how the crown "catalyzes" the guest's inversion. Two close-lying minima are characterized, as well as two distinct transition states connecting them, either via H3O+ inversion or rotation, with barriers of 1.0 and 4.6 kcal/mol, respectively, at the BLYP/II'//BLYP/6-31G level. Orbital interactions between lone pairs on ether O atoms and hydronium sigma(OH) antibonding orbitals are important factors for the directionality of the hydrogen bonds.     

Total synthesis of the microtubule stabilizing antitumor agent laulimalide and some nonnatural analogues: the power of Sharpless' asymmetric epoxidation

Three different routes are described for the synthesis of deoxylaulimalide (3), which is the immediate precursor of the marine sponge metabolite laulimalide (1). These routes mainly differ with respect to their ring closing step. Thus, route 1 uses a Still-Gennari olefination, route 2 a Yamaguchi lactonization, and route 3 an intramolecular allylsilane-aldehyde addition for establishing the macrocyclic structure. The unprotected deoxy derivative 3 was subjected to Sharpless' asymmetric epoxidation (SAE). With (R,R)-tartrate the 16,17-epoxide laulimalide (1) is formed selectively, whereas (S,S)-tartrate generates the 21,22-epoxide 142. This demonstrates the high reagent control involved in the SAE process, which in this case is used to achieve high stereo- and regioselectivity. Laulimalide and some derivatives thereof have been tested with respect to antitumor activity and compared to standard compounds paclitaxel and epothilone B.     

The C-Cl···π interactions inside supramolecular nanotubes of hexaethynylhexamethoxy[6]pericyclyne

Model complexes of hexaethynylhexamethoxy[6]pericyclyne (HM6P) molecules, with or without dichloromethane (DCM) guests, have been calculated at the M05-2X/6-311G(d,p) DFT level. Analysis of nonbonding interactions shows that the cohesion of the supramolecular tubular assemblies previously observed in the crystal state, relies mainly on C-H···O hydrogen bonds between axial ethynyl and equatorial methoxy substituents of stacked HM6P C(18) macrocycles in a cyclohexane-like chair conformation. The intrinsic binding energy of one HM6P molecule with two neighbors is calculated to be more than 40 kcal mol(-1). The inner channel of the stacks (of ca. 8 ? diameter) are suitable for hosting DCM molecules. Using the Atoms In Molecules (AIM) theory, the interaction between DCM molecules and surrounding triple bonds is analyzed in terms of σ-hole-directed C-Cl···π halogen bonds. A σ-hole-directed Cl···Cl interaction between two chlorine atoms of different DCM molecules is also evidenced.     

Is there a difference between leads and drugs? A historical perspective

To be considered for further development, lead structures should display the following properties: (1) simple chemical features, amenable for chemistry optimization; (2) membership to an established SAR series; (3) favorable patent situation; and (4) good absorption, distribution, metabolism, and excretion (ADME) properties. There are two distinct categories of leads: those that lack any therapeutic use (i.e., "pure" leads), and those that are marketed drugs themselves but have been altered to yield novel drugs. We have previously analyzed the design of leadlike combinatorial libraries starting from 18 lead and drug pairs of structures (S. J. Teague et al. Angew. Chem., Int. Ed. Engl. 1999, 38, 3743-3748). Here, we report results based on an extended dataset of 96 lead-drug pairs, of which 62 are lead structures that are not marketed as drugs, and 75 are drugs that are not presumably used as leads. We examined the following properties: MW (molecular weight), CMR (the calculated molecular refractivity), RNG (the number of rings), RTB (the number of rotatable bonds), the number of hydrogen bond donors (HDO) and acceptors (HAC), the calculated logarithm of the n-octanol/water partition (CLogP), the calculated logarithm of the distribution coefficient at pH 7.4 (LogD(74)), the Daylight-fingerprint druglike score (DFPS), and the property and pharmacophore features score (PPFS). The following differences were observed between the medians of drugs and leads: DeltaMW = 69; DeltaCMR = 1.8; DeltaRNG = DeltaHAC =1; DeltaRTB = 2; DeltaCLogP = 0.43; DeltaLogD(74) = 0.97; DeltaHDO = 0; DeltaDFPS = 0.15; DeltaPPFS = 0.12. Lead structures exhibit, on the average, less molecular complexity (less MW, less number of rings and rotatable bonds), are less hydrophobic (lower CLogP and LogD(74)), and less druglike (lower druglike scores). These findings indicate that the process of optimizing a lead into a drug results in more complex structures. This information should be used in the design of novel combinatorial libraries that are aimed at lead discovery.     

Molecular Dynamics,Phase Transition and Frequency‐Tuned Dielectric Switch of an Ionic Co‐Crystal

Dielectric switches that can be converted between high and low dielectric states by thermal stimuli have attracted much interest owing to their many potential applications. Currently one main drawback for practical application lies in the non‐tunability of their switch temperatures (T_S). We report here an ionic co‐crystal (Me₃NH)₄[Ni(NCS)₆] that contains a multiply rotatable Me₃NH⁺ ion and a solely rotatable one due to a more spacious supramolecular cage for the former one. This compound undergoes an isostructural order–disorder phase transition and it can function as a frequency‐tuned dielectric switch with highly adjustable T_S, which is further revealed by the variable‐temperature structure analyses and molecular dynamics simulations. In addition, the distinct arrangements and molecular dynamics of two coexisting Me₃NH⁺ ions confined in different lattice spaces as well as the notable offset effect on the promoting/hindering of dipolar reorientation after dielectric transition provide a rarely observed but fairly good model for understanding and modulating the dipole motion in crystalline environment.