Langevin Dipoles Model for ab Initio Calculations of Chemical Processes in Solution: Parametrization and Application to Hydration Free Energies of Neutral and Ionic Solutes and Conformational Analysis in Aqueous Solution

The Determination of Free Energy of Hydration of Water Ions from First Principles

We model the autoionization of water by determining the free energy of hydration of the major intermediate species of water ions. We represent the smallest ions─the hydroxide ion OH-, the hydronium ion H3O+, and the Zundel ion H5O2+─by bonded models and the more extended ionic structures by strong nonbonded interactions (e.g., the Eigen H9O4+ = H3O+ + 3(H2O) and the Stoyanov H13O6+ = H5O2+ + 4(H2O)). Our models are faithful to the precise QM energies and their components to within 1% or less. Using the calculated free energies and atomization energies, we compute the pKa of pure water from first principles as a consistency check and arrive at a value within 1.3 log units of the experimental one. From these calculations, we conclude that the hydronium ion, and its hydrated state, the Eigen cation, are the dominant species in the water autoionization process.

Solvation energies of the ferrous ion in water and in ammonia at various temperatures

CONTEXT

The solvation of metal ions is crucial to understanding relevant properties in physics, chemistry, or biology. Therefore, we present solvation enthalpies and solvation free energies of the ferrous ion in water and ammonia. Our results agree well with the experimental reports for the hydration free energy and hydration enthalpy. We obtained [Formula: see text] kJ mol[Formula: see text] for the hydration free energy and [Formula: see text] kJ mol[Formula: see text] for the hydration enthalpy of ferrous ion in water at room temperature. At ambient temperature, we obtained [Formula: see text] kJ mol[Formula: see text] as the [Formula: see text] ammoniation free energy and [Formula: see text] kJ mol[Formula: see text] for the ammoniation enthalpy. In addition, the free energy of solvation is deeply affected when the temperature increases. This pattern can be attributed to the rise of entropy when the temperature rises. Besides, the temperature does not affect the ammoniation enthalpies and the hydration enthalpy of the [Formula: see text] ion.

METHOD

All the geometry optimizations are performed at the MP2 methods associated with the 6-31++g(d,p) basis set of Pople. solvated phase structures of [Formula: see text] ion in water or in ammonia are performed using the PCM model. The [Formula: see text] program suite was used to perform all the calculations. The program TEMPO was also used to evaluate the temperature sensitivity of the different obtained geometries.

Ion Solvation Free Energy Calculation Based on Ab Initio Molecular Dynamics Using a Hybrid Solvent Model

Free energy calculation of small molecules or ion species in aqueous solvent is one of the most important problems in electrochemistry study. Although there are many previous approaches to calculate such free energies, they are based on ab initio methods and suffer from various limitations and approximations. In the current work, we developed a hybrid approach based on ab initio molecular dynamics (AIMD) simulations to calculate the ion solvation energy. In this approach, a small water cluster surrounding the central ion is used, and implicit solvent model is applied outside the water cluster. A dynamic potential well is used during AIMD to keep the water cluster together. Quasi-harmonic approximation is used to calculate the entropy contribution, while the total energy average is used to calculate the enthalpy term. The obtained solvation voltages of the bulk metal agree with experiments within 0.3 eV, and the simulation results for the solvation energies of gaseous ions are close to the experimental observations. Besides the free energies, radial pair distribution functions and coordination numbers of hydrated cations are also obtained. The remaining challenges of this method are also discussed.

Hydrolysis of Acetamide on Low-Index CeO2 Surfaces: Ceria as a Deamidation and General De-esterification Catalyst

Using DFT calculations and acetamide as the main example, we show that ceria is a potential catalyst for the hydrolysis of amide and similar bonds. The overall reaction is endergonic in the gas phase, yielding acetic acid and ammonia, but is slightly exergonic in the aqueous phase, which facilitates ionization of the products (CH₃COO^– and NH₄⁺). Neighboring Ce and O sites on the CeO₂(111), (110), and (100) facets are conducive to the formation of an activated metastable tetrahedral intermediate (TI) complex, followed by C–N bond scission. With van der Waals and solvation effects taken into account, the overall reaction energetics is found to be most favorable on the (111) facet as desorption of acetic acid is much more uphill energetically on (110) and (100). We further suggest that the Ce–O–Ce sites on ceria surfaces can activate X(=Y)–Z type bonds in amides, amidines, and carboxylate and phosphate esters, among many others that we term “generalized esters”. A Brønsted-Evans–Polanyi relationship is identified correlating the stability of the transition and final states of the X–Z generalized ester bond scission. A simple descriptor (ΣΔχ) based on the electronegativity of the atoms that constitute the bond (X, Y, Z) versus those of the catalytic site (O, Ce, Ce) captures the trend in the stability of the transition state of generalized ester bond scission and suggests a direction for modifying ceria for targeting specific organic substrates.

Simulating the binding of key organic functional groups to aqueous calcium carbonate species

The interaction of organic molecules with mineral systems is relevant to a wide variety of scientific problems both in the environment and minerals processing. In this study, the coordination of small organics that contain the two most relevant functional groups for biomineralisation of calcium carbonate, namely carboxylate and ammonium, with the corresponding mineral ions are examined in aqueous solution. Specifically, two force fields have been examined based on rigid-ion or polarisable models, with the latter being within the AMOEBA formalism. Here the parameters for the rigid-ion model are determined to target the accurate reproduction of the hydration structure and solvation thermodynamics, while both force fields are designed to be compatible with the corresponding recently published models for aqueous calcium carbonate. The application of these force fields to ion pairing in aqueous solution is studied in order to quantitatively determine the extent of association.

Calculating Entropies of Large Molecules in Aqueous Phase

Entropy benchmarking of different sized molecules in aqueous phase is carried out for known solvation models, where we compare geometry and solvation cavity packing parameters, which allows us to improve the accuracy of the obtained entropy values using empirical corrections. A comparison of solvation entropy models is conducted for a benchmarking set of 56 molecules, showing how an accurate description of cavitation entropy and its hindrance on other entropy values is important for large-sized solute molecules. Finally, we compare reaction free energies with entropies calculated using the most accurate solvation model considered, where we demonstrate a significant improvement in the accuracy relative to experimental values.

Chemoprotective Effects of Xanthohumol against the Carcinogenic Mycotoxin Aflatoxin B1

The present study addresses the chemoprotective effects of xanthohumol (XN), a prenylated flavonoid found in the female inflorescences (hops) of the plant Humulus lupulus L., against the carcinogenic food contaminant aflatoxin B1 (AFB1). The chemical reactions of XN and its derivatives (isoxanthohumol (IXN), 8-prenylnaringenin (8-PN), and 6-prenylnaringenin (6-PN)) with the AFB1 metabolite, aflatoxin B1 exo-8,9-epoxide (AFBO), were investigated in silico, by calculating activation free energies (ΔG^‡) at the Hartree–Fock level of theory in combination with the 6-311++G(d,p) basis set and two implicit solvation models. The chemoprotective effects of XN were investigated in vitro in the metabolically competent HepG2 cell line, analyzing its influence on AFB1-induced cytotoxicity using the MTS assay, genotoxicity using the comet and γH2AX assays, and cell cycle modulation using flow cytometry. Our results show that the ΔG^‡ required for the reactions of XN and its derivatives with AFBO are comparable to the ΔG^‡ required for the reaction of AFBO with guanine, indicating that XN, IXN, 8-PN, and 6-PN could act as scavengers of AFBO, preventing DNA adduct formation and DNA damage induction. This was also reflected in the results from the in vitro experiments, where a reduction in AFB1-induced cytotoxicity and DNA single-strand and double-strand breaks was observed in cells exposed to combinations of AFB1 and XN, highlighting the chemoprotective effects of this phytochemical.

Cyclization step of noradrenaline and adrenaline autoxidation: a quantum chemical study

Catecholamine autoxidation has been recognized as one of the potential trigger factors for catecholaminergic neuron loss characteristics of neurodegenerative diseases. The cyclization step with intramolecular Michael addition of catecholamine o-quinones has been shown to be the irreversible and rate limiting step of the autoxidation reaction across a broad pH range and has a complex pH dependence that has not yet been fully understood. Using quantum chemical calculations, we demonstrated that in the case of noradrenaline and adrenaline two catecholamine o-quinone species, one with an unprotonated and one with a protonated quinone group can participate in the cyclization reaction and that the mechanisms of these reactions significantly differ, emphasizing the importance of quinone group protonation states in the reaction mechanism. With a thorough exploration of the reaction kinetics, we further showed that at acidic pH the cyclization reaction rate is pH independent, while at alkaline pH the pH dependence is marked, explaining the experimentally observed complex pH dependence.

The quinone group protonation state determines the reaction mechanism of noradrenaline and adrenaline o-quinone cyclization.

Carcinogenesis of β-Propiolactone: A Computational Study

The discovery that β-propiolactone (BPL), once a commercially important chemical, causes various tumors in experimental animals has led to a significant decrease in its use. However, owing to its efficacy this possible human carcinogen remains to be utilized in vaccines for inactivation of viruses. The focus of the current study was to uncover the mechanisms of β-propiolactone reactions with both nucleobases and glutathione (GSH) through computer simulations based on quantum chemical methods. Our results, in accordance with in vitro studies, show that among all nucleobases guanine most readily forms adducts with BPL through S_N2 reaction mechanism. Acquired activation energies with incorporated solvent effects reveal that alkylation represents an energetically more favorable reaction than acylation for all nucleobases. Comparison of activation free energies of glutathione and guanine reactions with BPL suggest that glutathione may represent an efficient natural scavenger of BPL. Therefore, glutathione present in the organism may provide protection to the DNA and thus prevent BPL’s genotoxicity, mutagenicity, and possibly even carcinogenicity.

Protective Effects of [6]-Gingerol Against Chemical Carcinogens: Mechanistic Insights

[6]-Gingerol from ginger has received considerable attention as a potential cancer therapeutic agent because of its chemopreventive and chemotherapeutic effects, as well as its safety. In the current study, we examined [6]-gingerol as a natural scavenger of nine ultimate chemical carcinogens to which we are frequently exposed: glycidamide, styrene oxide, aflatoxin B1 exo-8,9-epoxide, β-propiolactone, ethylene oxide, propylene oxide, 2-cyanoethylene oxide, chloroethylene oxide, and vinyl carbamate epoxide. To evaluate [6]-gingerol efficacy, we expanded our research with the examination of glutathione-the strongest natural scavenger in human cells. The corresponding activation free energies were calculated using Hartree-Fock method with three flexible basis sets and two implicit solvation models. According to our results, [6]-gingerol proves to be an extremely effective scavenger of chemical carcinogens of the epoxy type. On the other hand, with the exception of aflatoxin B1 exo-8,9-epoxide, glutathione represents a relatively poor scavenger, whose efficacy could be augmented by [6]-gingerol. Moreover, our quantum mechanical study of the alkylation reactions of chemical carcinogens with [6]-gingerol and glutathione provide valuable insights in the reaction mechanisms and the geometries of the corresponding transition states. Therefore, we strongly believe that our research forms a solid basis for further computational, experimental and clinical studies of anticarcinogenic properties of [6]-gingerol as well as for the development of novel chemoprophylactic dietary supplements. Finally, the obtained results also point to the applicability of quantum chemical methods to studies of alkylation reactions related to chemical carcinogenesis.

Tannin Basic Building Blocks as Potential Scavengers of Chemical Carcinogens: A Computational Study

Tannins are natural compounds that have historically been used in the tanning of leather. In the scientific literature, one finds many reports of their possible beneficial health effects, although these are not always unequivocally confirmed. In order to gain a better insight into their proposed anticancer potential, we studied the scavenging capacity of the basic tannin building blocks against various chemical carcinogens of the epoxy type. The reactivity of gallic acid, ellagic acid, and epicathechin was examined using quantum mechanical calculations at the Hartree-Fock level of theory in conjunction with flexible basis sets and implicit solvation models. The monomeric tannin building blocks exhibited significant scavenging potential, with epicatechin presenting the best scavenger, thus encouraging and guiding future experimental studies of their anticarcinogenic properties.

Assessment of pKa Determination for Monocarboxylic Acids with an Accurate Theoretical Composite Method: G4CEP

For 22 monoprotic acids, pK_a values were calculated using the G4CEP composite theory. The solvation effect was included considering the continuous SMD solvation model; SMD and one explicit water molecule; and SMD, one water molecule, and linear correction with respect to the experimental pK_a values. The three tests provided mean absolute errors equal to 0.83, 0.51, and 0.30 pK_a units, respectively, indicating excellent performance of the G4CEP method. Comparison with density functional theory at the B3LYP and BMK levels showed that these results are quickly obtained but with a significant error. The best performance of the functionals was obtained from the combination of SMD, one explicit water molecule, linear regression correction, and basis set including diffuse functions. However, the dispersion of the results with DFT can lead to deviations of up to two pK_a units, whereas for G4CEP the largest deviations seldom exceed one pK_a unit.

Development of constant-pH simulation methods in implicit solvent and applications in biomolecular systems

pH is a critical parameter for biological and technological systems directly related with electrical charges. It can give rise to peculiar electrostatic phenomena, which also makes them more challenging. Due to the quantum nature of the process, involving the forming and breaking of chemical bonds, quantum methods should ideally by employed. Nevertheless, due to the very large number of ionizable sites, different macromolecular conformations, salt conditions, and all other charged species, the CPU time cost simply becomes prohibitive for computer simulations, making this a quite complex problem. Simplified methods based on Monte Carlo sampling have been devised and will be reviewed here, highlighting the updated state-of-the-art of this field, advantages, and limitations of different theoretical protocols for biomolecular systems (proteins and nucleic acids). Following a historical perspective, the discussion will be associated with the applications to protein interactions with other proteins, polyelectrolytes, and nanoparticles.

How fast monoamine oxidases decompose adrenaline? Kinetics of isoenzymes A and B evaluated by empirical valence bond simulation

This work scrutinizes kinetics of decomposition of adrenaline catalyzed by monoamine oxidase (MAO) A and B enzymes, a process controlling the levels of adrenaline in the central nervous system and other tissues. Experimental kinetic data for MAO A and B catalyzed decomposition of adrenaline are reported only in the form of the maximum reaction rate. Therefore, we estimated the experimental free energy barriers form the kinetic data of closely related systems using regression method, as was done in our previous study. By using multiscale simulation on the Empirical Valence Bond (EVB) level, we studied the chemical reactivity of the MAO A catalyzed decomposition of adrenaline and we obtained a value of activation free energy of 17.3 ± 0.4 kcal/mol. The corresponding value for MAO B is 15.7 ± 0.7 kcal/mol. Both values are in good agreement with the estimated experimental barriers of 16.6 and 16.0 kcal/mol for MAO A and MAO B, respectively. The fact that we reproduced the kinetic data and preferential catalytic effect of MAO B over MAO A gives additional support to the validity of the proposed hydride transfer mechanism. Furthermore, we demonstrate that adrenaline is preferably involved in the reaction in a neutral rather than in a protonated form due to considerably higher barriers computed for the protonated adrenaline substrate. The results are discussed in the context of chemical mechanism of MAO enzymes and possible applications of multiscale simulation to rationalize the effects of MAO activity on adrenaline level.

Averaged Condensed Phase Model for Simulating Molecules in Complex Environments

The need for configurational sampling dramatically increases the cost of combined quantum mechanics/molecular mechanics (QM/MM) simulations of chemical processes in solution. We developed an averaged condensed phase environment (ACPE) model that constructs an effective polarizable environment directly from explicitly sampled molecular dynamics configurations via the K-means++ algorithm and a mathematically rigorous translation of the molecular mechanics parameters. The model captures detailed heterogeneous features in the environment that may be difficult to describe using a conventional polarizable continuum model. Instead of performing repeated QM/MM calculations for each new configuration of the environment, the ACPE approach allows one to perform a single QM calculation on an averaged configuration. Here, we demonstrate the model by computing electronic excitation energies for several small molecules in solution. The ACPE model predicts the excitation energies in excellent agreement with conventional configurational averaging yet with orders of magnitude of reduction in the computational cost.

SM6: A Density Functional Theory Continuum Solvation Model for Calculating Aqueous Solvation Free Energies of Neutrals, Ions, and Solute-Water Clusters

A new charge model, called Charge Model 4 (CM4), and a new continuum solvent model, called Solvation Model 6 (SM6), are presented. Using a database of aqueous solvation free energies for 273 neutrals, 112 ions, and 31 ion-water clusters, parameter sets for the mPW0 hybrid density functional of Adamo and Barone (Adamo, C.; Barone, V. J. Chem. Phys. 1998, 108, 664-675) were optimized for use with the following four basis sets:  MIDI!6D, 6-31G(d), 6-31+G(d), and 6-31+G(d,p). SM6 separates the observable aqueous solvation free energy into two different components:  one arising from long-range bulk electrostatic effects and a second from short-range interactions between the solute and solvent molecules in the first solvation shell. This partition of the observable solvation free energy allows SM6 to effectively model a wide range of solutes. For the 273 neutral solutes in the test set, SM6 achieves an average error of ∼0.50 kcal/mol in the aqueous solvation free energies. For solutes, especially ions, that have highly concentrated regions of charge density, adding an explicit water molecule to the calculation significantly improves the performance of SM6 for predicting solvation free energies. The performance of SM6 was tested against several other continuum models, including SM5.43R and several different implementations of the Polarizable Continuum Model (PCM). For both neutral and ionic solutes, SM6 outperforms all of the models against which it was tested. Also, SM6 is the only model (except for one with an average error 3.4 times larger) that improves when an explicit solvent molecule is added to solutes with concentrated charge densities. Thus, in SM6, unlike the other continuum models tested here, adding one or more explicit solvent molecules to the calculation is an effective strategy for improving the prediction of the aqueous solvation free energies of solutes with strong local solute-solvent interactions. This is important, because local solute-solvent interactions are not specifically accounted for by bulk electrostatics, but modeling these interactions correctly is important for predicting the aqueous solvation free energies of certain solutes. Finally, SM6 retains its accuracy when used in conjunction with the B3LYP and B3PW91 functionals, and in fact the solvation parameters obtained with a given basis set may be used with any good density functional or fraction of Hartree-Fock exchange.

First Principles Calculations of Aqueous pKa Values for Organic and Inorganic Acids Using COSMO-RS Reveal an Inconsistency in the Slope of the pKa Scale

The COSMO-RS method, a combination of the quantum chemical dielectric continuum solvation model COSMO with a statistical thermodynamics treatment for more realistic solvation (RS) simulations, has been used for the direct prediction of pKa constants of a large variety of 64 organic and inorganic acids. A highly significant correlation of r(2) = 0.984 with a standard deviation of only 0.49 between the calculated values of the free energies of dissociation and the experimental pKa values was found, without any special adjustment of the method. Thus, we have a theoretical a priori prediction method for pKa, which has the regression constant and the slope as only adjusted parameters. Such a method can be of great value in many areas of physical chemistry, especially in pharmaceutical and agrochemical industry. To our surprise, the slope of pKa vs ΔGdiss is only 58% of the theoretically expected value of 1/RTln(10). A careful analysis with respect to different contributions as well as a comparison with the work of other authors excludes the possibility that the discrepancy is due to weaknesses of the calculation method. Hence, we must conclude that the experimental pKa scale depends differently on the free energy of dissociation than generally assumed.

Hydration in discrete water (II): from neutral to charged solutes

In our previous work, we introduced a solvation model based on discrete solvent representation and demonstrated its ability to estimate hydration free energies for neutral solutes. Here, we present modifications extending the applicability of the model to charged solutes. They include improvements in the representation of the first hydration shell and systematic treatment of long-range interactions. While sharing computational efficiency of implicit solvent models, our approach avoids some of their important limitations, both in the context of electrostatic and nonpolar hydration effects: it naturally captures hydration asymmetry of opposite charges, it relies on solute description by standard all atom force fields instead of utilizing specialized sets of atomic parameters, it predicts solvent distribution in space without the need to geometrically define solvent accessible surface. By combining discrete solvent representation in vicinity of a solute with continuum description of long-range interactions, the model addresses two distinct aspects of biomolecular hydration: complex, short-range effects arising due to molecular nature of aqueous solvent, and bulk contributions. We demonstrate that the model provides good agreement with experimental results for an extensive set of roughly 700 diverse compounds, including neutral and charged solutes with hydration free energies ranging from +3.4 to -536 kcal/mol.

Polarizable molecular interactions in condensed phase and their equivalent nonpolarizable models

Earlier, using phenomenological approach, we showed that in some cases polarizable models of condensed phase systems can be reduced to nonpolarizable equivalent models with scaled charges. Examples of such systems include ionic liquids, TIPnP-type models of water, protein force fields, and others, where interactions and dynamics of inherently polarizable species can be accurately described by nonpolarizable models. To describe electrostatic interactions, the effective charges of simple ionic liquids are obtained by scaling the actual charges of ions by a factor of 1/√(ε(el)), which is due to electronic polarization screening effect; the scaling factor of neutral species is more complicated. Here, using several theoretical models, we examine how exactly the scaling factors appear in theory, and how, and under what conditions, polarizable Hamiltonians are reduced to nonpolarizable ones. These models allow one to trace the origin of the scaling factors, determine their values, and obtain important insights on the nature of polarizable interactions in condensed matter systems.

Carcinogenesis of urethane: simulation versus experiment

The carcinogenesis of urethane (ethyl carbamate), a byproduct of fermentation that is consistently found in various food products, was investigated with a combination of kinetic experiments and quantum chemical calculations. The main objective of the study was to find ΔG(⧧), the activation free energy for the rate-limiting step of the SN2 reaction among the ultimate carcinogen of urethane, vinyl carbamate epoxide (VCE), and different nucleobases of the DNA. In the experimental part, the second-order reaction rate constants for the formation of the main 7-(2-oxoethyl)guanine adduct in aqueous solutions of deoxyguanosine and in DNA were determined. A series of ab initio, density functional theory (DFT), and semiempirical molecular orbital (MO) calculations was then performed to determine the activation barriers for the reaction between VCE and nucleobases methylguanine, methyladenine, and methylcytosine. Effects of hydration were incorporated with the use of the solvent reaction field method of Tomasi and co-workers and the Langevine dipoles model of Florian and Warshel. The computational results for the main adduct were found to be in good agreement with the experiment, thus presenting strong evidence for the validity of the proposed SN2 mechanism. This allowed us to predict the activation barriers of reactions leading to side products for which kinetic experiments have not yet been performed. Our calculations have shown that the main 7-(2-oxoethyl)deoxyguanosine adduct indeed forms preferentially because the emergence of other adducts either proceeds across a significantly higher activation barrier or the geometry of the reaction requires the Watson-Crick pairs of the DNA to be broken. The computational study also considered the questions of stereoselectivity, the ease of the elimination of the leaving group, and the relative contributions of the two possible reaction paths for the formation of the 1,N(2)-ethenodeoxyguanosine adduct.

Recent advances in QM/MM free energy calculations using reference potentials

Background

Recent years have seen enormous progress in the development of methods for modeling (bio)molecular systems. This has allowed for the simulation of ever larger and more complex systems. However, as such complexity increases, the requirements needed for these models to be accurate and physically meaningful become more and more difficult to fulfill. The use of simplified models to describe complex biological systems has long been shown to be an effective way to overcome some of the limitations associated with this computational cost in a rational way.

Scope of review

Hybrid QM/MM approaches have rapidly become one of the most popular computational tools for studying chemical reactivity in biomolecular systems. However, the high cost involved in performing high-level QM calculations has limited the applicability of these approaches when calculating free energies of chemical processes. In this review, we present some of the advances in using reference potentials and mean field approximations to accelerate high-level QM/MM calculations. We present illustrative applications of these approaches and discuss challenges and future perspectives for the field.

Major conclusions

The use of physically-based simplifications has shown to effectively reduce the cost of high-level QM/MM calculations. In particular, lower-level reference potentials enable one to reduce the cost of expensive free energy calculations, thus expanding the scope of problems that can be addressed.

General significance

As was already demonstrated 40 years ago, the usage of simplified models still allows one to obtain cutting edge results with substantially reduced computational cost. This article is part of a Special Issue entitled Recent developments of molecular dynamics.

•

We present some of the advances to accelerate high-level QM/MM calculations.

•

Quantitative limitations of low-level methods can be overcome by these approaches.

•

Reference potentials make free energy simulations feasible for large systems.

•

Automated fitting reduces the need of expensive sampling of high-level approaches.

•

Application of reference potentials can be extended to a wide range of processes.

Quantum mechanical analysis of nonenzymatic nucleotidyl transfer reactions: kinetic and thermodynamic effects of β-γ bridging groups of dNTP substrates

Rate (k) and equilibrium (K) constants for the reaction of tetrahydrofuranol with a series of Mg²⁺ complexes of methyl triphosphate analogues, CH₃O-P(O₂)-O-P(O₂)-X-PO₃^4–, X = O, CH₂, CHCH₃, C(CH₃)₂, CFCH₃, CHF, CHCl, CHBr, CFCl, CF₂, CCl₂, and CBr₂, forming phosphate diester and pyrophosphate or bisphosphonate in aqueous solution were evaluated by B3LYP/TZVP//HF/6-31G* quantum chemical calculations and Langevin dipoles and polarized continuum solvation models. The calculated log k and log K values were found to depend linearly on the experimental pK_a4 of the conjugate acid of the corresponding pyrophosphate or bisphosphonate leaving group. The calculated slopes of these Brønsted linear free energy relationships were β_lg = −0.89 and β_eq = −0.93, respectively. The studied compounds also followed the linear relationship Δlog k = 0.8Δlog K, which became less steep, Δlog k = 0.6Δlog K, after the range of studied compounds was extended to include analogues that were doubly protonated on γ-phosphate, CH₃O-P(O₂)-O-P(O₂)-X-PO₃H₂^2–. The scissile P_α–O_lg bond length in studied methyl triphosphate analogues slightly increases with decreasing pK_a of the leaving group; concomitantly, the CH₃OP_α(O₂) moiety becomes more positive. These structural effects indicate that substituents with low pK_a can facilitate both P_α–O_lg bond breaking and the P_α–O_nuc bond forming process, thus explaining the large negative β_lg calculated for the transition state geometry that has significantly longer P_α–O_nuc distance than the P_α–O_lg distance.

Hydration free energies of cyanide and hydroxide ions from molecular dynamics simulations with accurate force fields

The evaluation of hydration free energies is a sensitive test to assess force fields used in atomistic simulations. We showed recently that the vibrational relaxation times, 1D- and 2D-infrared spectroscopies for CN(-) in water can be quantitatively described from molecular dynamics (MD) simulations with multipolar force fields and slightly enlarged van der Waals radii for the C- and N-atoms. To validate such an approach, the present work investigates the solvation free energy of cyanide in water using MD simulations with accurate multipolar electrostatics. It is found that larger van der Waals radii are indeed necessary to obtain results close to the experimental values when a multipolar force field is used. For CN(-), the van der Waals ranges refined in our previous work yield hydration free energy between -72.0 and -77.2 kcal mol(-1), which is in excellent agreement with the experimental data. In addition to the cyanide ion, we also study the hydroxide ion to show that the method used here is readily applicable to similar systems. Hydration free energies are found to sensitively depend on the intermolecular interactions, while bonded interactions are less important, as expected. We also investigate in the present work the possibility of applying the multipolar force field in scoring trajectories generated using computationally inexpensive methods, which should be useful in broader parametrization studies with reduced computational resources, as scoring is much faster than the generation of the trajectories.

Mechanisms of enzyme-catalyzed reduction of two carcinogenic nitro-aromatics, 3-nitrobenzanthrone and aristolochic acid I: Experimental and theoretical approaches

This review summarizes the results found in studies investigating the enzymatic activation of two genotoxic nitro-aromatics, an environmental pollutant and carcinogen 3-nitrobenzanthrone (3-NBA) and a natural plant nephrotoxin and carcinogen aristolochic acid I (AAI), to reactive species forming covalent DNA adducts. Experimental and theoretical approaches determined the reasons why human

NAD(P)H

quinone oxidoreductase (NQO1) and cytochromes P450 (CYP) 1A1 and 1A2 have the potential to reductively activate both nitro-aromatics. The results also contributed to the elucidation of the molecular mechanisms of these reactions. The contribution of conjugation enzymes such as N,O-acetyltransferases (NATs) and sulfotransferases (SULTs) to the activation of 3-NBA and AAI was also examined. The results indicated differences in the abilities of 3-NBA and AAI metabolites to be further activated by these conjugation enzymes. The formation of DNA adducts generated by both carcinogens during their reductive activation by the NOQ1 and CYP1A1/2 enzymes was investigated with pure enzymes, enzymes present in subcellular cytosolic and microsomal fractions, selective inhibitors, and animal models (including knock-out and humanized animals). For the theoretical approaches, flexible in silico docking methods as well as ab initio calculations were employed. The results summarized in this review demonstrate that a combination of experimental and theoretical approaches is a useful tool to study the enzyme-mediated reaction mechanisms of 3-NBA and AAI reduction.

Modified Poisson–Boltzmann equations for characterizing biomolecular solvation

Methods for computing electrostatic interactions often account implicitly for the solvent, due to the much smaller number of degrees of freedom involved. In the Poisson–Boltzmann (PB) approach the electrostatic potential is obtained by solving the Poisson–Boltzmann equation (PBE), where the solvent region is modeled as a homogeneous medium with a high dielectric constant. PB however is not exempt of problems. It does not take into account for example the sizes of the ions in the atmosphere surrounding the solute, nor does it take into account the inhomogeneous dielectric response of water due to the presence of a highly charged surface. In this paper we review two major modifications of PB that circumvent these problems, namely the size-modified PB (SMPB) equation and the Dipolar Poisson–Boltzmann Langevin (DPBL) model. In SMPB, steric effects between ions are accounted for with a lattice gas model. In DPBL, the solvent region is no longer modeled as a homogeneous dielectric media but rather as an assembly of self-orienting interacting dipoles of variable density. This model results in a dielectric profile that transits smoothly from the solute to the solvent region as well as in a variable solvent density that depends on the charges of the solute. We show successful applications of the DPBL formalism to computing the solvation free energies of isolated ions in water. Further developments of more accurately modified PB models are discussed.

Electrostatic-Consistent Coarse-Grained Potentials for Molecular Simulations of Proteins

We present a new generation of coarse-grained (CG) potentials that account for a simplified electrostatic description of soluble proteins. The treatment of permanent electrostatic dipoles of the backbone and polar side-chains allows to simulate proteins, preserving an excellent structural and dynamic agreement with respective reference structures and all-atom molecular dynamics simulations. Moreover, multiprotein complexes can be well described maintaining their molecular interfaces thanks to the ability of this scheme to better describe the actual electrostatics at a CG level of resolution. An efficient and robust heuristic algorithm based on particle swarm optimization is used for the derivation of CG parameters via a force-matching procedure. The ability of this protocol to deal with high dimensional search spaces suggests that the extension of this optimization procedure to larger data sets may lead to the generation of a fully transferable CG force field. At the present stage, these electrostatic-consistent CG potentials are easily and efficiently parametrized, show a good degree of transferability, and can be used to simulate soluble proteins or, more interestingly, large macromolecular assemblies for which long all-atom simulations may not be easily affordable.

Computational chemical analysis of unconjugated bilirubin anions and insights into pKa values clarification

The pKa, the negative logarithm of the acid dissociation equilibrium constant, of the carboxylic acid groups of unconjugated bilirubin in water is a discussed issue because there are quite different experimental values reported. Using quantum mechanical calculations we have studied the conformational behavior of unconjugated bilirubin species (in gas phase and in solution modeled implicitly and explicitly) to provide evidence that may clarify pKa values because of its pathophysiological relevance. Our results show that rotation of carboxylate group, which is not restricted, settles it in a suitable place to establish stronger interactions that stabilizes the monoanion and the dianion to be properly solvated, demonstrating that the rationalization used to justify the high pKa values of unconjugated bilirubin is inappropriate. Furthermore, low unconjugated bilirubin (UCB) pKa values were estimated from a linear regression analysis.

A multiscale coarse-grained polarizable solvent model for handling long tail bulk electrostatics

A multiscale coarse-grained approach able to handle efficiently the solvation of microscopic solutes in extended chemical environment is described. That approach is able to compute readily and efficiently very long-range solute/solvent electrostatic microscopic interactions, up to the 1-μm scale, by considering a reduced amount of computational resources. All the required parameters are assigned to reproduce available data concerning the solvation of single ions. Such a strategy makes it possible to reproduce with good accuracy the solvation properties concerning simple ion pairs in solution (in particular, the asymptotic behavior of the ion pair potentials of mean force). This new method represents an extension of the polarizable pseudoparticle solvent model, which has been recently improved to account for the main features of hydrophobic effects in liquid water (Masella et al., J. Comput. Chem. 2011, 32, 2664). This multiscale approach is well suited to be used for computing the impact of charge changes in free energy computations, in terms of both accuracy and efficiency.

Biomolecular electrostatics and solvation: a computational perspective

An understanding of molecular interactions is essential for insight into biological systems at the molecular scale. Among the various components of molecular interactions, electrostatics are of special importance because of their long-range nature and their influence on polar or charged molecules, including water, aqueous ions, proteins, nucleic acids, carbohydrates, and membrane lipids. In particular, robust models of electrostatic interactions are essential for understanding the solvation properties of biomolecules and the effects of solvation upon biomolecular folding, binding, enzyme catalysis, and dynamics. Electrostatics, therefore, are of central importance to understanding biomolecular structure and modeling interactions within and among biological molecules. This review discusses the solvation of biomolecules with a computational biophysics view toward describing the phenomenon. While our main focus lies on the computational aspect of the models, we provide an overview of the basic elements of biomolecular solvation (e.g. solvent structure, polarization, ion binding, and non-polar behavior) in order to provide a background to understand the different types of solvation models.

Spectroscopic and DFT studies of second-sphere variants of the type 1 copper site in azurin: covalent and nonlocal electrostatic contributions to reduction potentials

The reduction potentials (E(0)) of type 1 (T1) or blue copper (BC) sites in proteins and enzymes with identical first coordination spheres around the redox active copper ion can vary by ~400 mV. Here, we use a combination of low-temperature electronic absorption and magnetic circular dichroism, electron paramagnetic resonance, resonance Raman, and S K-edge X-ray absorption spectroscopies to investigate a series of second-sphere variants--F114P, N47S, and F114N in Pseudomonas aeruginosa azurin--which modulate hydrogen bonding to and protein-derived dipoles nearby the Cu-S(Cys) bond. Density functional theory calculations correlated to the experimental data allow for the fractionation of the contributions to tuning E(0) into covalent and nonlocal electrostatic components. These are found to be significant, comparable in magnitude, and additive for active H-bonds, while passive H-bonds are mostly nonlocal electrostatic in nature. For dipoles, these terms can be additive to or oppose one another. This study provides a methodology for uncoupling covalency from nonlocal electrostatics, which, when coupled to X-ray crystallographic data, distinguishes specific local interactions from more long-range protein/active interactions, while affording further insight into the second-sphere mechanisms available to the protein to tune the E(0) of electron-transfer sites in biology.