Skip to main page content
U.S. flag

An official website of the United States government

Dot gov

The .gov means it’s official.
Federal government websites often end in .gov or .mil. Before sharing sensitive information, make sure you’re on a federal government site.

Https

The site is secure.
The https:// ensures that you are connecting to the official website and that any information you provide is encrypted and transmitted securely.

Access keys NCBI Homepage MyNCBI Homepage Main Content Main Navigation
Review
. 2009 Feb 5;113(5):1253-72.
doi: 10.1021/jp8071712.

Progress in ab initio QM/MM free-energy simulations of electrostatic energies in proteins: accelerated QM/MM studies of pKa, redox reactions and solvation free energies

Affiliations
Review

Progress in ab initio QM/MM free-energy simulations of electrostatic energies in proteins: accelerated QM/MM studies of pKa, redox reactions and solvation free energies

Shina C L Kamerlin et al. J Phys Chem B. .

Abstract

Hybrid quantum mechanical/molecular mechanical (QM/MM) approaches have been used to provide a general scheme for chemical reactions in proteins. However, such approaches still present a major challenge to computational chemists, not only because of the need for very large computer time in order to evaluate the QM energy but also because of the need for proper computational sampling. This review focuses on the sampling issue in QM/MM evaluations of electrostatic energies in proteins. We chose this example since electrostatic energies play a major role in controlling the function of proteins and are key to the structure-function correlation of biological molecules. Thus, the correct treatment of electrostatics is essential for the accurate simulation of biological systems. Although we will be presenting different types of QM/MM calculations of electrostatic energies (and related properties) here, our focus will be on pKa calculations. This reflects the fact that pKa's of ionizable groups in proteins provide one of the most direct benchmarks for the accuracy of electrostatic models of macromolecules. While pKa calculations by semimacroscopic models have given reasonable results in many cases, existing attempts to perform pKa calculations using QM/MM-FEP have led to discrepancies between calculated and experimental values. In this work, we accelerate our QM/MM calculations using an updated mean charge distribution and a classical reference potential. We examine both a surface residue (Asp3) of the bovine pancreatic trypsin inhibitor and a residue buried in a hydrophobic pocket (Lys102) of the T4-lysozyme mutant. We demonstrate that, by using this approach, we are able to reproduce the relevant side chain pKa's with an accuracy of 3 kcal/mol. This is well within the 7 kcal/mol energy difference observed in studies of enzymatic catalysis, and is thus sufficient accuracy to determine the main contributions to the catalytic energies of enzymes. We also provide an overall perspective of the potential of QM/MM calculations in general evaluations of electrostatic free energies, pointing out that our approach should provide a very powerful and accurate tool to predict the electrostatics of not only solution but also enzymatic reactions, as well as the solvation free energies of even larger systems, such as nucleic acid bases incorporated into DNA.

PubMed Disclaimer

Figures

Figure 1
Figure 1
Energy scheme for evaluating the solvation free energy. Here, QISM0 is the charge distribution obtained by a given implicit solvation model, ΔGsol(QISM0) is evaluated by the classical FEP-AC approach, and all other terms are evaluated using a hybrid QM/MM approach.
Figure 2
Figure 2
Model for evaluating the average solute (S) charges for (a) a molecule in solution and (b) a protein sidechain. All atoms in Region I are represented explicitly, while Region II atoms are represented by two charges. In case (a), Region I only comprises explicit water molecules while in case (b) it also comprises all electroneutral groups within a pre-defined cut-off of the relevant protein sidechain. “P” designates the protein.
Figure 3
Figure 3
A graphical description of the thermodynamic cycle for estimating the energetics of acid dissociation in a protein. The ionisation process has been represented relative to the energetics and solvation of the A and AH species for the corresponding dissociation in water. Here, ΔGp and ΔGw represent the free energy of ionising an acid in a protein and water respectively, ΔGsolvpwandΔGsolvwp represent the difference in solvation energy of the indicated species in protein (p) and water (w) respectively, R is the ideal gas constant, and AH and A represent the ionised an neutral forms of the acid respectively.
Figure 4
Figure 4
The position of the Asp3 sidechain on the surface of the bovine pancreatic trypsin inhibitor (BPTI). The protein is shown in violet, and Asp3 is coloured by atom type. Also highlighted are all water molecules and electroneutral groups (shown in magenta) within 5Å of Asp3.
Figure 5
Figure 5
The contribution to the total free energy of solvation from the LRA terms of Eq. 6. (a) shows the free contribution for acetate in water, and (b) shows contribution for the Asp3 sidechain of the bovine pancreatic trypsin inhibitor (BPTI). In both cases, the neutral form is shown in blue and the ionised form in magenta.
Figure 6
Figure 6
The position of the Lys102 sidechain in the M102K lysozyme mutant. The protein is shown in violet, and Lys is coloured by atom type. Also highlighted are all water molecules and electroneutral groups (shown in magenta) within 5Å of Lys102.
Figure 7
Figure 7
The contribution to the total free energy of solvation from the LRA terms of Eq. 6. (a) shows the free contribution for methylamine in water, and (b) shows contribution for the Lys102 sidechain of the M102K mutated T4-lysozyme. In both cases, the neutral form is shown in blue and the ionised form in magenta.
Figure 8
Figure 8
Average (a) solute-solvent interaction energy and (b) solute polarization energy, along a 10 ps simulation of HCOO solvated by water (using the same solvation model as that presented in Fig. 2b). The fact that the final sum of <Eint> and <Epol> is independent of m cannot be realised from this figure, but rather from Fig. 9 of Ref..
Figure 9
Figure 9
Free energy surface for the hydrolysis of the methyl phosphate dianion (MeOPO3 2−). Here, RS denotes the reactant state, PS the product state, and TS1 and TS2 associative and dissociative transition states respectively.
Figure 10
Figure 10
A comparison of reaction profiles obtained by the COSMO solvation model (red) and the QM/MM-FEP approach presented in this work (blue) for hydroxide attack on the 4-nitro substituted methyl phenyl phosphate diester shown to the right, as a function of distance between the phosphorus atom being attacked and the oxygen atom of the incoming nucleophile (P-Onuc distance). Distances are given in Å, and energies (relative to the diester and hydroxide at infinite separation) are given in kcal/mol.

Similar articles

Cited by

References

    1. Cui Q, Elstner M, Kaxiras E, Frauenheim T, Karplus M. J. Phys. Chem. B. 2001;105:569–585.
    1. Monard G, Merz KM. Acc. Chem. Res. 1999;32:904–911.
    1. Shurki A, Warshel A. Adv. Protein. Chem. 2003;66:249–312. - PubMed
    1. Warshel A, Levitt M. J. Mol. Biol. 1976;103:227–249. - PubMed
    1. Gao J. Acc. Chem. Res. 1996;29:298–305.

Publication types