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.

Neural Network Corrections to Intermolecular Interaction Terms of a Molecular Force Field Capture Nuclear Quantum Effects in Calculations of Liquid Thermodynamic Properties.

We incorporate nuclear quantum effects (NQE) in condensed matter simulations by introducing short-range neural network (NN) corrections to the ab initio fitted molecular force field ARROW. Force field NN corrections are fitted to average interaction energies and forces of molecular dimers, which are simulated using the Path Integral Molecular Dynamics (PIMD) technique with restrained centroid positions. The NN-corrected force field allows reproduction of the NQE for computed liquid water and methane properties such as density, radial distribution function (RDF), heat of evaporation (HVAP), and solvation free energy. Accounting for NQE through molecular force field corrections circumvents the need for explicit computationally expensive PIMD simulations in accurate calculations of the properties of chemical and biological systems. The accuracy and locality of pairwise NN NQE corrections indicate that this approach could be applicable to complex heterogeneous systems, such as proteins.

Combining Force Fields and Neural Networks for an Accurate Representation of Bonded Interactions.

We present a formalism of a neural network encoding bonded interactions in molecules. This intramolecular encoding is consistent with the models of intermolecular interactions previously designed by this group. Variants of the encoding fed into a corresponding neural network may be used to economically improve the representation of torsional degrees of freedom in any force field. We test the accuracy of the reproduction of the ab initio potential energy surface on a set of conformations of two dipeptides, methyl-capped ALA and ASP, in several scenarios. The encoding, either alone or in conjunction with an analytical potential, improves agreement with ab initio energies that are on par with those of other neural network-based potentials. Using the encoding and neural nets in tandem with an analytical model places the agreements firmly within "chemical accuracy" of ±0.5 kcal/mol.

Combining Force Fields and Neural Networks for an Accurate Representation of Chemically Diverse Molecular Interactions.

A key goal of molecular modeling is the accurate reproduction of the true quantum mechanical potential energy of arbitrary molecular ensembles with a tractable classical approximation. The challenges are that analytical expressions found in general purpose force fields struggle to faithfully represent the intermolecular quantum potential energy surface at close distances and in strong interaction regimes; that the more accurate neural network approximations do not capture crucial physics concepts, e.g., nonadditive inductive contributions and application of electric fields; and that the ultra-accurate narrowly targeted models have difficulty generalizing to the entire chemical space. We therefore designed a hybrid wide-coverage intermolecular interaction model consisting of an analytically polarizable force field combined with a short-range neural network correction for the total intermolecular interaction energy. Here, we describe the methodology and apply the model to accurately determine the properties of water, the free energy of solvation of neutral and charged molecules, and the binding free energy of ligands to proteins. The correction is subtyped for distinct chemical species to match the underlying force field, to segment and reduce the amount of quantum training data, and to increase accuracy and computational speed. For the systems considered, the hybrid ab initio parametrized Hamiltonian reproduces the two-body dimer quantum mechanics (QM) energies to within 0.03 kcal/mol and the nonadditive many-molecule contributions to within 2%. Simulations of molecular systems using this interaction model run at speeds of several nanoseconds per day.

Protein-Ligand Binding Free-Energy Calculations with ARROW─A Purely First-Principles Parameterized Polarizable Force Field.

Protein-ligand binding free-energy calculations using molecular dynamics (MD) simulations have emerged as a powerful tool for in silico drug design. Here, we present results obtained with the ARROW force field (FF)─a multipolar polarizable and physics-based model with all parameters fitted entirely to high-level ab initio quantum mechanical (QM) calculations. ARROW has already proven its ability to determine solvation free energy of arbitrary neutral compounds with unprecedented accuracy. The ARROW FF parameterization is now extended to include coverage of all amino acids including charged groups, allowing molecular simulations of a series of protein-ligand systems and prediction of their relative binding free energies. We ensure adequate sampling by applying a novel technique that is based on coupling the Hamiltonian Replica exchange (HREX) with a conformation reservoir generated via potential softening and nonequilibrium MD. ARROW provides predictions with near chemical accuracy (mean absolute error of ∼0.5 kcal/mol) for two of the three protein systems studied here (MCL1 and Thrombin). The third protein system (CDK2) reveals the difficulty in accurately describing dimer interaction energies involving polar and charged species. Overall, for all of the three protein systems studied here, ARROW FF predicts relative binding free energies of ligands with a similar accuracy level as leading nonpolarizable force fields.

Accurate determination of solvation free energies of neutral organic compounds from first principles.

The main goal of molecular simulation is to accurately predict experimental observables of molecular systems. Another long-standing goal is to devise models for arbitrary neutral organic molecules with little or no reliance on experimental data. While separately these goals have been met to various degrees, for an arbitrary system of molecules they have not been achieved simultaneously. For biophysical ensembles that exist at room temperature and pressure, and where the entropic contributions are on par with interaction strengths, it is the free energies that are both most important and most difficult to predict. We compute the free energies of solvation for a diverse set of neutral organic compounds using a polarizable force field fitted entirely to ab initio calculations. The mean absolute errors (MAE) of hydration, cyclohexane solvation, and corresponding partition coefficients are 0.2 kcal/mol, 0.3 kcal/mol and 0.22 log units, i.e. within chemical accuracy. The model (ARROW FF) is multipolar, polarizable, and its accompanying simulation stack includes nuclear quantum effects (NQE). The simulation tools' computational efficiency is on a par with current state-of-the-art packages. The construction of a wide-coverage molecular modelling toolset from first principles, together with its excellent predictive ability in the liquid phase is a major advance in biomolecular simulation.

Harnessing the Substrate Promiscuity of Dioxygenase AsqJ and Developing Efficient Chemoenzymatic Synthesis for Quinolones.

Nature has developed complexity-generating reactions within natural product biosynthetic pathways. However, direct utilization of these pathways to prepare compound libraries remains challenging due to limited substrate scopes, involvement of multiple-step reactions, and moderate robustness of these sophisticated enzymatic transformations. Synthetic chemistry, on the other hand, offers an alternative approach to prepare natural product analogs. However, owing to complex and diverse functional groups appended on the targeted molecules, dedicated design and development of synthetic strategies are typically required. Herein, by leveraging the power of chemo-enzymatic synthesis, we report an approach to bridge the gap between biological and synthetic strategies in the preparation of quinolone alkaloid analogs. Leading by in silico analysis, the predicted substrate analogs were chemically synthesized. The AsqJ-catalyzed asymmetric epoxidation of these substrate analogues was followed by an Lewis Acid-triggered ring contraction to complete the viridicatin formation. We evaluated the robustness of this method in gram-scale reactions. Lastly, through chemoenzymatic cascades, a library of quinolone alkaloids is effectively prepared.

Epoxidation Catalyzed by the Nonheme Iron(II)- and 2-Oxoglutarate-Dependent Oxygenase, AsqJ: Mechanistic Elucidation of Oxygen Atom Transfer by a Ferryl Intermediate.

Mechanisms of enzymatic epoxidation via oxygen atom transfer (OAT) to an olefin moiety is mainly derived from the studies on thiolate-heme containing epoxidases, such as cytochrome P450 epoxidases. The molecular basis of epoxidation catalyzed by nonheme-iron enzymes is much less explored. Herein, we present a detailed study on epoxidation catalyzed by the nonheme iron(II)- and 2-oxoglutarate-dependent (Fe/2OG) oxygenase, AsqJ. The native substrate and analogues with different para substituents ranging from electron-donating groups (e.g., methoxy) to electron-withdrawing groups (e.g., trifluoromethyl) were used to probe the mechanism. The results derived from transient-state enzyme kinetics, Mössbauer spectroscopy, reaction product analysis, X-ray crystallography, density functional theory calculations, and molecular dynamic simulations collectively revealed the following mechanistic insights: (1) The rapid O2 addition to the AsqJ Fe(II) center occurs with the iron-bound 2OG adopting an online-binding mode in which the C1 carboxylate group of 2OG is trans to the proximal histidine (His134) of the 2-His-1-carboxylate facial triad, instead of assuming the offline-binding mode with the C1 carboxylate group trans to the distal histidine (His211); (2) The decay rate constant of the ferryl intermediate is not strongly affected by the nature of the para substituents of the substrate during the OAT step, a reactivity behavior that is drastically different from nonheme Fe(IV)-oxo synthetic model complexes; (3) The OAT step most likely proceeds through a stepwise process with the initial formation of a C(benzylic)-O bond to generate an Fe-alkoxide species, which is observed in the AsqJ crystal structure. The subsequent C3-O bond formation completes the epoxide installation.

Phospholipid-Based Prodrugs for Colon-Targeted Drug Delivery: Experimental Study and In-Silico Simulations.

In ulcerative colitis (UC), the inflammation is localized in the colon, and one of the successful strategies for colon-targeting drug delivery is the prodrug approach. In this work, we present a novel phospholipid (PL)-based prodrug approach, as a tool for colonic drug targeting in UC. We aim to use the phospholipase A2 (PLA2), an enzyme that is overexpressed in the inflamed colonic tissues of UC patients, as the PL-prodrug activating enzyme, to accomplish the liberation of the parent drug from the prodrug complex at the specific diseased tissue(s). Different linker lengths between the PL and the drug moiety can dictate the rate of activation by PLA2, and subsequently determine the amount of free drugs at the site of action. The feasibility of this approach was studied with newly synthesized PL-Fmoc (fluorenylmethyloxycarbonyl) conjugates, using Fmoc as a model compound for testing our hypothesis. In vitro incubation with bee venom PLA2 demonstrated that a 7-carbon linker between the PL and Fmoc has higher activation rate than a 5-carbon linker. 4-fold higher colonic expression of PLA2 was demonstrated in colonic mucosa of colitis-induced rats when compared to healthy animals, validating our hypothesis of a colitis-targeting prodrug approach. Next, a novel molecular dynamics (MD) simulation was developed for PL-based prodrugs containing clinically relevant drugs. PL-methotrexate conjugate with 6-carbon linker showed the highest extent of PLA2-mediated activation, whereas shorter linkers were activated to a lower extent. In conclusion, this work demonstrates that for carefully designed PL-drug conjugates, PLA2 overexpression in inflamed colonic tissues can be used as prodrug-activating enzyme and drug targeting strategy, including insights into the activation mechanisms in a PLA2 binding site.

On the importance of accounting for nuclear quantum effects in ab initio calibrated force fields in biological simulations.

In many important processes in chemistry, physics, and biology the nuclear degrees of freedom cannot be described using the laws of classical mechanics. At the same time, the vast majority of molecular simulations that employ wide-coverage force fields treat atomic motion classically. In light of the increasing desire for and accelerated development of quantum mechanics (QM)-parameterized interaction models, we reexamine whether the classical treatment is sufficient for a simple but crucial chemical species: alkanes. We show that when using an interaction model or force field in excellent agreement with the "gold standard" QM data, even very basic simulated properties of liquid alkanes, such as densities and heats of vaporization, deviate significantly from experimental values. Inclusion of nuclear quantum effects via techniques that treat nuclear degrees of freedom using the laws of classical mechanics brings the simulated properties much closer to reality.

Insights into the Desaturation of Cyclopeptin and its C3 Epimer Catalyzed by a non-Heme Iron Enzyme: Structural Characterization and Mechanism Elucidation.

AsqJ, an iron(II)- and 2-oxoglutarate-dependent enzyme found in viridicatin-type alkaloid biosynthetic pathways, catalyzes sequential desaturation and epoxidation to produce cyclopenins. Crystal structures of AsqJ bound to cyclopeptin and its C3 epimer are reported. Meanwhile, a detailed mechanistic study was carried out to decipher the desaturation mechanism. These findings suggest that a pathway involving hydrogen atom abstraction at the C10 position of the substrate by a short-lived FeIV -oxo species and the subsequent formation of a carbocation or a hydroxylated intermediate is preferred during AsqJ-catalyzed desaturation.

Crystal Structure of Cleaved Serp-1, a Myxomavirus-Derived Immune Modulating Serpin: Structural Design of Serpin Reactive Center Loop Peptides with Improved Therapeutic Function.

The Myxomavirus-derived protein Serp-1 has potent anti-inflammatory activity in models of vasculitis, lupus, viral sepsis, and transplant. Serp-1 has also been tested successfully in a Phase IIa clinical trial in unstable angina, representing a "first-in-class" therapeutic. Recently, peptides derived from the reactive center loop (RCL) have been developed as stand-alone therapeutics for reducing vasculitis and improving survival in MHV68-infected mice. However, both Serp-1 and the RCL peptides lose activity in MHV68-infected mice after antibiotic suppression of intestinal microbiota. Here, we utilize a structure-guided approach to design and test a series of next-generation RCL peptides with improved therapeutic potential that is not reduced when the peptides are combined with antibiotic treatments. The crystal structure of cleaved Serp-1 was determined to 2.5 Å resolution and reveals a classical serpin structure with potential for serpin-derived RCL peptides to bind and inhibit mammalian serpins, plasminogen activator inhibitor 1 (PAI-1), anti-thrombin III (ATIII), and α-1 antitrypsin (A1AT), and target proteases. Using in silico modeling of the Serp-1 RCL peptide, S-7, we designed several modified RCL peptides that were predicted to have stronger interactions with human serpins because of the larger number of stabilizing hydrogen bonds. Two of these peptides (MPS7-8 and -9) displayed extended activity, improving survival where activity was previously lost in antibiotic-treated MHV68-infected mice (P < 0.0001). Mass spectrometry and kinetic assays suggest interaction of the peptides with ATIII, A1AT, and target proteases in mouse and human plasma. In summary, we present the next step toward the development of a promising new class of anti-inflammatory serpin-based therapeutics.

Computational modeling and in-vitro/in-silico correlation of phospholipid-based prodrugs for targeted drug delivery in inflammatory bowel disease.

Targeting drugs to the inflamed intestinal tissue(s) represents a major advancement in the treatment of inflammatory bowel disease (IBD). In this work we present a powerful in-silico modeling approach to guide the molecular design of novel prodrugs targeting the enzyme PLA2, which is overexpressed in the inflamed tissues of IBD patients. The prodrug consists of the drug moiety bound to the sn-2 position of phospholipid (PL) through a carbonic linker, aiming to allow PLA2 to release the free drug. The linker length dictates the affinity of the PL-drug conjugate to PLA2, and the optimal linker will enable maximal PLA2-mediated activation. Thermodynamic integration and Weighted Histogram Analysis Method (WHAM)/Umbrella Sampling method were used to compute the changes in PLA2 transition state binding free energy of the prodrug molecule (∆∆Gtr) associated with decreasing/increasing linker length. The simulations revealed that 6-carbons linker is the optimal one, whereas shorter or longer linkers resulted in decreased PLA2-mediated activation. These in-silico results were shown to be in excellent correlation with experimental in-vitro data. Overall, this modern computational approach enables optimization of the molecular design of novel prodrugs, which may allow targeting the free drug specifically to the diseased intestinal tissue of IBD patients.

Prediction of cyclohexane-water distribution coefficient for SAMPL5 drug-like compounds with the QMPFF3 and ARROW polarizable force fields.

We present the performance of blind predictions of water-cyclohexane distribution coefficients for 53 drug-like compounds in the SAMPL5 challenge by three methods currently in use within our group. Two of them utilize QMPFF3 and ARROW, polarizable force-fields of varying complexity, and the third uses the General Amber Force-Field (GAFF). The polarizable FF's are implemented in an in-house MD package, Arbalest. We find that when we had time to parametrize the functional groups with care (batch 0), the polarizable force-fields outperformed the non-polarizable one. Conversely, on the full set of 53 compounds, GAFF performed better than both QMPFF3 and ARROW. We also describe the torsion-restrain method we used to improve sampling of molecular conformational space and thus the overall accuracy of prediction. The SAMPL5 challenge highlighted several drawbacks of our force-fields, such as our significant systematic over-estimation of hydrophobic interactions, specifically for alkanes and aromatic rings.

Phospholipid-Based Prodrugs for Drug Targeting in Inflammatory Bowel Disease: Computational Optimization and In-Vitro Correlation.

In inflammatory bowel disease (IBD) patients, the enzyme phospholipase A2 (PLA2) is overexpressed in the inflamed intestinal tissue, and hence may be exploited as a prodrug-activating enzyme allowing drug targeting to the site(s) of gut inflammation. The purpose of this work was to develop powerful modern computational approaches, to allow optimized a-priori design of phospholipid (PL) based prodrugs for IBD drug targeting. We performed simulations that predict the activation of PL-drug conjugates by PLA2 with both human and bee venom PLA2. The calculated results correlated well with in-vitro experimental data. In conclusion, a-priori drug design using a computational approach complements and extends experimentally derived data, and may improve resource utilization and speed drug development.

Modeling Electronic Polarizability Changes in the Course of a Magnesium Ion Water Ligand Exchange Process.

This paper introduces explicit dependence of atomic polarizabilities on intermolecular interactions within the framework of a polarizable force field AMOEBA. Polarizable models used in biomolecular simulations often poorly describe molecular electrostatic induction in condensed phase, in part, due to neglect of a strong dependency of molecular electronic polarizability on intermolecular interactions at short distances. Our variable polarizability model parameters are derived from quantum chemical calculations of small clusters of atoms and molecules, and can be applied in simulations in condensed phase without additional scaling factors. The variable polarizability model is applied to simulate a ligand exchange reaction for a Mg(2+) ion solvated in water. Explicit dependence of water polarizability on a distance between a water oxygen and Mg(2+) is derived from in vacuum MP2 calculations of Mg(2+)-water dimer. The simulations yield a consistent description of the energetics of the Mg(2+)-water clusters of different size. Simulations also reproduce thermodynamics of ion solvation as well as kinetics of a water ligand exchange reaction. In contrast, simulations that used the additive force field or that used the constant polarizability models were not able to consistently and quantitatively describe the properties of the solvated Mg(2+) ion.

Donor/acceptor coupling shortcuts in electron transfer within ruthenium-modified derivatives of cytochrome b(562).

Quantitative theoretical studies of long-range electron transfer are still rare, and reliable computational methods to analyze these reactions are still being developed. We re-examined electron transfer reactions in ruthenium-modified cytochrome b562 derivatives focusing on accurate calculation of statistical average of electron transfer rates that are dominated by a small fraction of accessible protein conformations. We performed a series of ab initio calculations of donor/acceptor interactions over protein fragments sampled from long molecular dynamic trajectories and compared computed electron transfer rates to available experimental data. Our approach takes into account cofactor electronic structure and effects of solvation on the donor-acceptor interactions. It allows predicting absolute values of electron transfer rates in contrast to other computational methodologies that give only qualitative results. Our calculations reproduced with a good accuracy experimental electron transfer rates. We also found that electron transfer in some of the cytochrome b562 derivatives is dominated by "shortcut" conformations, where donor/acceptor interactions are mediated by nonbonded interactions of Ru ligands with protein surface groups. Several derivatives adopt long-lived conformations with the Ru complex interacting with negatively charged protein residues that are characterized by shorter Ru-Fe distances and higher ET rates. We argue that quantitative theoretical analysis is essential for detailed understanding of protein electron transfer and mechanisms of biological redox reactions.

A glance at the reactivity of osma(II)cycles [Os(C-N)x(bpy)3-x](m+) (x=0-3) Covering a 1.8 V potential range toward peroxidase through Monte Carlo Simulations ((-)C-N=o-2-phenylpyridinato, bpy=2,2'-bipyridine).

Three cyclometalated and one coordination compounds [Os(C-N)x(bpy)3-x](m) (x/m=0/2+ (4); 1/1+ (3); 2/1+ (2); 3/0 (1); (-)C-N=2-phenylpyridinato, bpy=2,2'-bipyridine) with drastically different reduction potentials have been used for analyzing the second-order rate constants for one-electron, metal-based osmium(II) to osmium(III) oxidation of the complexes by compound I (k2) and compound II (k3) of horseradish peroxidase. Previously unknown k2 and k3 have been determined by digital simulation of cyclic voltammograms measured in phosphate buffer of pH7.6 and 21 ± 1°C. Osmium(II) species derived from osmium(III) complexes 1 and 2 were generated electrochemically in situ. Under the conditions used the reduction potentials for the Os(III/II) feature equal -0.90, -0.095, 0.23 and 0.85V versus NHE (normal hydrogen electrode) for 1-4, respectively. The rate constants k2 equal ~5 × 10(7), 6 × 10(8), 2 × 10(6) and 1 × 10(5)M(-1)s(-1) and the rate constants k3 equal ~9 × 10(6), 4× 10(7), 1 ×10(6) and 1 × 10(5)M(-1)s(-1) for complexes 1-4, respectively. Both rate constants k2 and k3 first increase with increasing the reaction driving force on going from 4 to 2 but then both decrease on going to complex 1 though the reaction driving force is the highest in this case. The system described has been explored theoretically using docking Monte Carlo simulations.

Control of cytochrome c redox reactivity through off-pathway modifications in the protein hydrogen-bonding network.

Measurements of photoinduced Fe(2+)-to-Ru(3+) electron transfer (ET), supported by theoretical analysis, demonstrate that mutations off the dominant ET pathways can strongly influence the redox reactivity of cytochrome c. The effects arise from the change in the protein dynamics mediated by the intraprotein hydrogen-bonding network.

pH-triggered conformational switching of the diphtheria toxin T-domain: the roles of N-terminal histidines.

pH-induced conformational switching is essential for functioning of diphtheria toxin, which undergoes a membrane insertion/translocation transition triggered by endosomal acidification as a key step of cellular entry. In order to establish the sequence of molecular rearrangements and side-chain protonation accompanying the formation of the membrane-competent state of the toxin's translocation (T) domain, we have developed and applied an integrated approach that combines multiple techniques of computational chemistry [e.g., long-microsecond-range, all-atom molecular dynamics (MD) simulations; continuum electrostatics calculations; and thermodynamic integration (TI)] with several experimental techniques of fluorescence spectroscopy. TI calculations indicate that protonation of H257 causes the greatest destabilization of the native structure (6.9 kcal/mol), which is consistent with our early mutagenesis results. Extensive equilibrium MD simulations with a combined length of over 8 µs demonstrate that histidine protonation, while not accompanied by the loss of structural compactness of the T-domain, nevertheless results in substantial molecular rearrangements characterized by the partial loss of secondary structure due to unfolding of helices TH1 and TH2 and the loss of close contact between the C- and N-terminal segments. The structural changes accompanying the formation of the membrane-competent state ensure an easier exposure of the internal hydrophobic hairpin formed by helices TH8 and TH9, in preparation for its subsequent transmembrane insertion.