Models of F.H contacts relevant to the binding of fluoroaromatic inhibitors to carbonic anhydrase II

[formula: see text] Complexes formed between fluorobenzene and N-methylformamide or benzene have been used as models of the interaction of fluoroaromatic drugs with carbonic anhydrase II. These structures have been investigated via ab initio and density functional methods, including HF, B3LYP, and MP2 procedures. The results of the calculations are consistent with the hypothesis, suggested originally by experimental X-ray crystal structures of the drug-receptor complexes, that favorable fluorine-hydrogen interactions affect binding affinity.     

Linear free energy relationships implicate three modes of binding for fluoroaromatic inhibitors to a mutant of carbonic anhydrase II

[figure: see text] Linear free energy relationships between binding affinity and hydrophobicity for a library of fluoroaromatic inhibitors of F131V carbonic anhydrase II (CA) implicate three modes of interaction. X-ray crystal structures suggest that F131 interacts with fluoroaromatic inhibitors, while P202, on the opposite side of the active site cleft, serves as the site of the hydrophobic contact in the case of the F131V mutant. 2-Fluorinated compounds bind more tightly, perhaps due to the field effect of the nearby fluorine on the acidity of the amide proton.     

Chaperoned alchemical free energy simulations: a general method for QM, MM, and QM/MM potentials

A general method for alchemical free energy simulations using QM, MM, and QM/MM potential is developed by introducing "chaperones" to restrain the structures, particularly near the end points. A calculation of the free energy difference between two triazole tautomers in aqueous solution is used to illustrate the method.     

pKa analysis for the zinc-bound water in human carbonic anhydrase II: Benchmark for "multiscale" QM/MM simulations and mechanistic implications

To quantitatively explore the applicability of the generalized solvent boundary potential (GSBP) based QM/MM approach as a "multiscale" framework for studying chemical reactions in biomolecules, the structural and energetic properties of the Human Carbonic Anhydrase II (CAII) are analyzed and compared to those from periodic boundary condition (PBC) simulations and available experimental data. Although the atomic fluctuations from GSBP based simulations are consistently lower compared to those from PBC simulations or crystallographic data, the fluctuations and internal coordinate distributions for residues in the proximity of the active site as well as diffusion constants of active-site water molecules are fairly well described by GSBP simulations. The pKa of the zinc-bound water, calculated with a SCC-DFTB/MM-GSBP based thermodynamic integration approach, agrees well with experiments for the wild type CAII. For the E106Q mutant, however, a 9 pKa unit downward shift relative to the wild type is found in contrast with previous experiments that found little change. This dramatic discrepancy signals a possible change in the mechanism for the interconversion between CO2/HCO3- in the E106Q mutant, which may be similar to the bicarbonate mediated mechanism proposed for the Co2+ substituted CAII (J. Am. Chem. Soc. 2001, 123, 5861).1 The study highlights pKa analyses as a valuable approach for quantitatively validating the computational model for complex biomolecules as well as for revealing energetic properties intimately related to the chemical process of interest.     

The pattern of fluorine substitution affects binding affinity in a small library of fluoroaromatic inhibitors for carbonic anhydrase

[formula: see text] A library of fluoroaromatic inhibitors of carbonic anhydrase has been found to bind in a manner dependent on both hydrophobicity and the pattern of substitution of the fluoroaromatic ring. All of the compounds in the library bind to the protein with Kd < 3 nM. We have inferred two distinct binding modes from our data, which suggest two types of interactions that should be considered when designing fluorinated drugs.     

A QM/MM study of the binding of RAPTA ligands to cathepsin B

We have carried out quantum mechanical (QM) and QM/MM (combined QM and molecular mechanics) calculations, as well as molecular dynamics (MD) simulations to study the binding of a series of six RAPTA (Ru(II)-arene-1,3,5-triaza-7-phosphatricyclo-[3.3.1.1] decane) complexes with different arene substituents to cathepsin B. The recently developed QM/MM-PBSA approach (QM/MM combined with Poisson–Boltzmann solvent-accessible surface area solvation) has been used to estimate binding affinities. The QM calculations reproduce the antitumour activities of the complexes with a correlation coefficient (r ²) of 0.35–0.86 after a conformational search. The QM/MM-PBSA method gave a better correlation (r ² = 0.59) when the protein was fixed to the crystal structure, but more reasonable ligand structures and absolute binding energies were obtained if the protein was allowed to relax, indicating that the ligands are strained when the protein is kept fixed. In addition, the best correlation (r ² = 0.80) was obtained when only the QM energies were used, which suggests that the MM and continuum solvation energies are not accurate enough to predict the binding of a charged metal complex to a charged protein. Taking into account the protein flexibility by means of MD simulations slightly improves the correlation (r ² = 0.91), but the absolute energies are still too large and the results are sensitive to the details in the calculations, illustrating that it is hard to obtain stable predictions when full flexible protein is included in the calculations.     

A new QM/MM method oriented to the study of ionic liquids

The interest on room temperature ionic liquids has grown in the last decades because of their use as all‐purpose solvent and their low environmental impact. In the present work, a new theoretical procedure is developed to study pure ionic liquids within the framework of the quantum mechanics/molecular mechanics method. Each type of ion (cation or anion) is considered as an independent entity quantum mechanically described that follows a differentiated path in the liquid. The method permits, through an iterative procedure, the full coupling between the polarized charge distribution of the ions and the liquid structure around them. The procedure has been tested with 1‐ethyl‐3‐methylimidazolium tetrafluoroborate. It was found that, similar to non‐polar liquids and as a consequence of the low value of the reaction field, the cation and anion charge distributions are hardly polarized by the rest of molecules in the liquid. Their structure is characterized by an alternance between anion and cation shells as evidenced by the coincidence of the first maximum of the anion–anion and cation–cation radial distribution functions with the first minimum of the anion‐cation. Some degree of stacking between the cations is also found. © 2015 Wiley Periodicals, Inc.     

Unraveling the origin of regioselectivity in rhodium diphosphine catalyzed hydroformylation. A DFT QM/MM study

The origin of regioselectivity in rhodium diphosphine catalyzed hydroformilation was investigated by means of hybrid QM/MM calculations using the IMOMM method. The roles of the diphosphine bite angle and of the nonbonding interactions were analyzed in detail by considering rhodium systems containing xantphos-type ligands, for which a correlation between the natural bite angle and regioselectivity has been recently reported. From the pentacoordinated equatorial--equatorial HRh(CO)(alkene)(diphosphine) key intermediate, eight possible reaction paths were defined and characterized through their respective transition states (TS). We succeeded in reproducing the experimentally observed trends for the studied diphosphines. By performing additional calculations on model systems, in which the steric effects induced by the phenyl substituents of xantphos ligands were canceled, we were able to separate, identify, and evaluate the different contributions to regioselectivity. These additional calculations showed that regioselectivity is governed by the nonbonding interactions between the diphenylphosphino substituents and the substrate, whereas the effects directly associated to the bite angle, what we call orbital effects, seem to have a smaller influence.     

Molecular dynamics study of a new rigidified 18-crown-6 derivative using a QM / MM method

We present results of the first molecular dynamics study of a new rigidified 18-crown-6 derivative, recently synthesized by Li and Still. We refer to this new crown as s18c6, which contains a central 18-crown-6 core locked into a D_3d configuration by the addition of six exocyclic saturated hydrocarbon rings, two of which also contain a single exocyclic heteroatom (X), peripheral to the central cation binding ring. X can be either a hemiacetal (X = oxygen) or hemithioacetal (X = sulfur). The enhanced rigidity and the exocyclic heteroatoms are thought to be involved in the novel cation binding properties of s18c6 in aqueous solutions; it is more ionophoric than 18c6, and it is selective for Na⁺ over K⁺ when X = oxygen but is selective for K⁺ over Na⁺ when X = sulfur. In the present communication, we report results for the X = oxygen isomer of s18c6. We use a hybrid quantum mechanical/molecular mechanical (QM/MM) approach to study the solvation properties of the uncomplexed crown in H₂O. The QM/MM method described here is based on our previous QM/MM study of K⁺/18c6 which employs the semiempirical AM1 method for s18c6 and the SPC/E model for H₂O. Both AM1 and HF/6–31 + G* optimized geometries give similar results for the gas-phase structure of s18c6. The ab initio results show that s186 has a slightly smaller cavity than 18c6 (5.681 vs. 5.802 Å, respectively) based on the average transannular oxygen-oxygen distances. A 100-ps equilibrium molecular dynamics simulation shows that the s18c6 core remains in a rigid D_3d configuration. S18c6 avidly binds two solvent waters, one on either side of the macrocycle plane. The simulation average binding energy per bound water is - 14.8 kcal/mol, which is greater than our previously reported result of - 12.3 kcal/mol for H₂O bound by 18c6. Both of the s18c6-bound waters maintain multiple hydrogen bonds to the oxygens of the macrocycle, and there is evidence for involvement of the exocyclic oxygens in this binding. © 1996 John Wiley & Sons, Inc.     

A critical evaluation of different QM/MM frontier treatments with SCC-DFTB as the QM method

The performance of different link atom based frontier treatments in QM/MM simulations was evaluated critically with SCC-DFTB as the QM method. In addition to the analysis of gas-phase molecules as in previous studies, an important element of the present work is that chemical reactions in realistic enzyme systems were also examined. The schemes tested include all options available in the program CHARMM for SCC-DFTB/MM simulation, which treat electrostatic interactions due to the MM atoms close to the QM/MM boundary in different ways. In addition, a new approach, the divided frontier charge (DIV), has been implemented in which the partial charge associated with the frontier MM atom ("link host") is evenly distributed to the other MM atoms in the same group. The performance of these schemes was evaluated based on properties including proton affinities, deprotonation energies, dipole moments, and energetics of proton transfer reactions. Similar to previous work, it was found that calculated proton affinities and deprotonation energies of alcohols, carbonic acids, amino acids, and model DNA bases are very sensitive to the link atom scheme; the commonly used single link atom approach often gives error on the order of 15 to 20 kcal/mol. Other schemes give better and, on average, mutually comparable results. For proton transfer reactions, encouragingly, both activation barriers and reaction energies are fairly insensitive (within a typical range of 2-4 kcal/mol) to the link atom scheme due to error cancellation, and this was observed for both gas-phase and enzyme systems. Therefore, the effect of using different link atom schemes in QM/MM simulations is rather small for chemical reactions that conserve the total charge. Although the current study used an approximate DFT method as the QM level, the observed trends are expected to be applicable to QM/MM methods with use of other QM approaches. This observation does not mean to encourage QM/MM simulations without careful benchmark in the study of specific systems, rather it emphasizes that other technical details, such as the treatment of long-range electrostatics, tend to play a more important role and need to be handled carefully.     

Theoretical studies on substrate binding to the active site of carbonic anhydrase

The role of some amino acids and metal ions in the catalytic activity of carbonic anhydrase (carbonate dehydratase EC 4.2.1.1) has been investigated. The additional stabilization of the transition state complex was used as the qualitative measure of the effect of molecular surrounding on CO₂ hydration reaction calculated within the approximate CNDO /2 approach. The effect of the molecular environment has been simulated by inclusion into the SCF LCAO MO Hamiltonian, a term representing the interaction with a set of point charges and atomic dipoles centered on experimentally determined position of Zn²⁺ ion and all atoms of histidine 94, 96, 119 and threonine 199, 200 (or histidine 200 in the case of carbonic anhydrase B). The possible molecular mechanism of CO₂ hydration inside the active site has been also discussed.     

QM/MM study of the taxadiene synthase mechanism

Combined quantum mechanics/molecular mechanics (QM/MM) calculations were used to investigate the reaction mechanism of taxadiene synthase (TXS). TXS catalyzes the cyclization of geranylgeranyl diphosphate (GGPP) to taxadiene (T) and four minor cyclic products. All these products originate from the deprotonation of carbocation intermediates. The reaction profiles for the conversion of GGPP to T as well as to minor products were calculated for different configurations of relevant TXS carbocation complexes. The QM region was treated at the M06-2X/TZVP level, while the CHARMM27 force field was used to describe the MM region. The QM/MM calculations suggest a reaction pathway for the conversion of GGPP to T, which slightly differs from previous proposals regarding the number of reaction steps and the conformation of the carbocations. The QM/MM results also indicate that the formation of minor products via water-assisted deprotonation of the carbocations is highly exothermic, by about −7 to −23 kcal/mol. Curiously, however, the computed barriers and reaction energies indicate that the formation of some of the minor products is more facile than the formation of T. Thus, the present QM/MM calculations provide detailed insights into possible reaction pathways and into the origin of the promiscuity of TXS, but they do not reproduce the product distribution observed experimentally. © 2019 Wiley Periodicals, Inc.     

Flexible-boundary QM/MM calculations: II. Partial charge transfer across the QM/MM boundary that passes through a covalent bond

Recently, based on the principle of electronic chemical potential equalization and the principle of charge conservation, we proposed a flexible-boundary scheme that allows both partial charge transfer and self-consistent polarization between the quantum mechanical (QM) and molecular mechanical (MM) subsystems in QM/MM calculations; the scheme was applied to study the atomic charges in selected ion–solvent complexes. In the present contribution, we further extend the flexible-boundary treatment to handle the QM/MM boundary passing through covalent bonds. We find that the flexible-boundary redistributed charge and dipole schemes yield reasonable agreement with full-QM calculations for a number of molecular ions and amino acids with charged side chains. Using the full-QM results as reference, the mean unsigned deviations are computed to be 0.06 e for atomic partial charges of the QM atoms, 0.11 e for the amounts of charge transfer between the QM and MM subsystems, and 0.016 Å for the lengths of the covalent bonds that directly connect the QM and MM subsystems. The results indicate the importance of accounting for partial charge transfer across the QM/MM boundary when the QM subsystems are charged.     

A QM/MM study on the catalytic mechanism of pyruvate decarboxylase

Pyruvate decarboxylase (PDC) is a typical thiamin diphosphate (ThDP)-dependent enzyme with widespread applications in industry. Though studies regarding the reaction mechanism of PDC have been reported, they are mainly focused on the formation of ThDP ylide and some elementary steps in the catalytic cycle, studies about the whole catalytic cycle of PDC are still not completed. In these previous studies, a major controversy is whether the key active residues (Glu473, Glu50′, Asp27′, His113′, His114′) are protonated or ionized during the reaction. To explore the catalytic mechanism and the role of key residues in the active site, three whole-enzyme models were considered, and the combined QM/MM calculations on the nonoxidative decarboxylation of pyruvate to acetaldehyde catalyzed by PDC were performed. According to our computational results, the fundamental reaction pathways, the complete energy profiles of the whole catalytic cycle, and the specific role of key residues in the common steps were obtained. It is also found that the same residue with different protonation states will lead to different reaction pathways and energy profiles. The mechanism derived from the model in which the residues (Glu473, Glu50′, Asp27′, His113′, His114′) are in their protonated states is most consistent with experimental observations. Therefore, extreme care must be taken when assigning the protonation states in the mechanism study. Because the experimental determination of protonation state is currently difficult, the combined QM/MM method provides an indirect means for determining the active-site protonation state.     

Theoretical modeling of the heme group with a hybrid QM/MM method

The quality of the results obtained in calculations with the hybrid QM/MM method IMOMM on systems where the heme group is partitioned in QM and MM regions is evaluated through the performance of calculations on the 4‐coordinate [Fe(P)] (P = porphyrin), the 5‐coordinate [Fe(P)(1−(Me)Im)] (Im = imidazole) and the 6‐coordinate [Fe(P)(1−(Me)Im)(O₂)] systems. The results are compared with those obtained from much more expensive pure quantum mechanics calculations on model systems. Three different properties are analyzed—namely, the optimized geometries, the binding energies of the axial ligands to the heme group, and the energy cost of the biochemically relevant out‐of‐plane displacement of the iron atom. Agreement is especially good in the case of optimized geometries and energy cost of out‐of‐plane displacements, with larger discrepancies in the case of binding energies. © 2000 John Wiley & Sons, Inc. J Comput Chem 21: 282–294, 2000     

A QM/MM simulation study of transamination reaction at the active site of aspartate aminotransferase: Free energy landscape and proton transfer pathways

Transaminase is a key enzyme for amino acid metabolism, which reversibly catalyzes the transamination reaction with the help of PLP (pyridoxal 5^' -phosphate) as its cofactor. Here we have investigated the mechanism and free energy landscape of the transamination reaction involving the aspartate transaminase (AspTase) enzyme and aspartate-PLP (Asp-PLP) complex using QM/MM simulation and metadynamics methods. The reaction is found to follow a stepwise mechanism where the active site residue Lys258 acts as a base to shuttle a proton from α -carbon (CA) to imine carbon (C4A) of the PLP-Asp Schiff base. In the first step, the Lys258 abstracts the CA proton of the substrate leading to the formation of a carbanionic intermediate which is followed by the reprotonation of the Asp-PLP Schiff base at C4A atom by Lys258. It is found that the free energy barrier for the proton abstraction by Lys258 and that for the reprotonation are 17.85 and 3.57 kcal/mol, respectively. The carbanionic intermediate is 7.14 kcal/mol higher in energy than the reactant. Hence, the first step acts as the rate limiting step. The present calculations also show that the Lys258 residue undergoes a conformational change after the first step of transamination reaction and becomes proximal to C4A atom of the Asp-PLP Schiff base to favor the second step. The active site residues Tyr70* and Gly38 anchor the Lys258 in proper position and orientation during the first step of the reaction and stabilize the positive charge over Lys258 generated at the intermediate step.     

The reaction mechanism of hydroxyethylphosphonate dioxygenase: a QM/MM study

By employing ab initio quantum mechanical/molecular mechanical (QM/MM) and molecular dynamics (MD) simulations, we have provided further evidence against the previously proposed hydroperoxylation or hydroxylation mechanism of hydroxyethylphosphonate dioxygenase (HEPD). HEPD employs an interesting catalytic cycle based on concatenated bifurcations. The first bifurcation is based on the abstraction of hydrogen atoms from the substrate, which leads to a distal or proximal hydroperoxo species (Fe-OOH or Fe-(OH)O). The second and the third bifurcations refer to the carbon-carbon bond cleavage reaction. And this is achieved through a tridentate intermediate, or employing a proton-shuttle assisted mechanism, in which the residue Glu(176) or the Fe(IV)=O group serves as a general base. The reaction directions seem to be tunable and show significant environment dependence. This mechanism can provide a comprehensive interpretation for the seemingly contradicting experimental evidences and provide insight into the development of biochemistry and material sciences.