Author: kararanaghan

In this recent publication included in the ‘Women in Computational Chemistry’ special issue of the Journal of Chemical Information and Modeling, we investigate the dependence of predicted reaction energetics on the choice of density functional and QM region size used in QM/MM calculations. DFT methods potentially offer a good combination of accuracy and computational cost but suffer from some well-known limitations and are not systematically improvable. The Claisen rearrangement of chorismate to prephenate, a simple unimolecular reaction catalysed by the enzyme chorismate mutase, was used as a model system to investigate how these choices for QM/MM calculations impact the predicted energetics of reaction.

There are many different density functionals available and the best density functional to choose for any application is not obvious. Commonly used density functionals (e.g. B3LYP, BH&HLYP, MO6-2X) predict barrier heights for the reaction in the enzyme that differ by as much 13 kcal/mol. When different density functionals give different results for the same system, which result should be preferred? Ab initio methods are potentially highly accurate and are systematically improvable but often restricted by their computational cost. Projector-based embedding allows the incorporation of ab initio methods into QM/MM calculations at reasonable computational cost as a small number of atoms can be selected for treatment at the highest levels. Embedding SCS-MP2 in DFT in these QM/MM calculations reduced the spread in predicted barrier height from 13 kcal/mol to 0.3 kcal/mol, essentially eliminating the dependence on the density functional, as shown below (top: standard DFT, bottom: SCS-MP2-in-DFT).

 

The optimum size of the QM region in QM/MM calculations is a matter of much debate in the literature. The effect of the size of the QM region was tested in solution (by adding additional water molecules to the QM region) and in the enzyme (by adding additional arginine side chains) showing that the difference in barrier heights predicted by common density functionals for any size of the QM region is significantly larger than the change in barrier caused by increasing the size of the QM region.

 

Ranaghan, K. E.,  Shchepanovska, D., Bennie, S. J.,  Lawan, N., Macrae, S. J., Zurek, J., Manby, F. R. and Mulholland, A. J.
J. Chem. Inf. Model., 2019
DOI: 10.1021/acs.jcim.8b00940

Constant pH Molecular Dynamics (CpHMD) is a simulation technique gaining in popularity owing to its ability to titrate ionisable residues and predict pK_a’s using standard molecular mechanics force fields. State of the art CpHMD simulations involve running pH replica exchange MD (pH-REMD) simulations using an explicit solvent model to accurately sample the protein conformational and protonation state spaces. It has been successfully applied to a variety of natural proteins to estimate the pK_a of key residues and provide mechanistic insights into pH-dependant phenomena.

In a recent paper published in Nature Communications, BCompB members Eric Lang and Adrian Mulholland pushed the boundaries of CpHMD simulations further by investigating the pH dependency of de novo α-helical barrels designed with a ring of 6 or 7 of strongly interacting glutamate residues whose side chains point toward the inside of the central, otherwise hydrophobic, lumen.

The CpHMD simulations were used to investigate the pH-dependant stability of the hexameric and heptameric α-helical barrels (named CC-Type2-LL-L17E and CC-Type2-IL-Sg-L17E respectively) observed experimentally by circular dichroism (CD). For each barrel, two 200 ns pH-REMD simulations in TIP3P water, were carried out with all ionisable residues allowed to titrate. A total of 16 replicas was used to cover a pH range between 3.0 and 10.5, with one replica per 0.5 pH unit, representing an aggregate sampling time of 6.4 µs for each structure.

Because the glutamate residues (Glu-17) are in close proximity and interact strongly, the authors anticipated that the individual microscopic pK_a of each glutamate residue made little sense and instead considered the Glu-17 rings as single polyprotic species composed of 6 or 7 titrable groups. Using this approach, a macroscopic stepwise pK_a can be calculated for the (Glu-17)₆ and (Glu-17)₇ species, giving results in excellent agreement with the experimental data. These calculations revealed that up to two negative charges can be accommodated inside the barrels but increasing the number of charged Glu dramatically decreases the stability of the barrels. The simulations also revealed that the two buried negative charges are accommodated due to the solvation of the rings of Glu-17 residues and the presence of a Na⁺ counter ion in the middle of the rings. Increasing the number of deprotonated Glu-17 above two leads to a reorientation of the Glu side chains to access bulk solvent. This, in turn, leads to disruption of the secondary structure and ultimately to the opening of the barrels.

The ring of Glu-17 residues in the central lumen of a hexameric α-helical barrel

 

The paper, which primarily presents the experimental and design work lead by Guto Rhys and Derek Woolfson, also provides an example of how the Isambard software, developed in the Wolfson group, was used by BCompB member Christopher Wood to rationalise the diversity of structures formed by the peptides studied in this work.

Growing resistance to antibiotics is a major threat to human health and there is an urgent need to develop new antibiotics or find new ways to keep existing antibiotics effective.  β-lactamases are a class of enzymes that have evolved to break down penicillin-type antibiotics, leading to resistance. One method to avoid this is to administer a β-lactamase inhibitor as a resistance blocker. Clavulanic acid is one such example and is often combined with amoxicillin to form the drug commonly known as co-amoxiclav. In a recent paper published in Biochemistry, van der Kamp and co-workers have used combined quantum mechanics / molecular mechanics (QM/MM) methods to investigate how clavulanic acid can be broken down by different β-lactamases. TEM-1 and KPC-2 are clinically important enzymes that exhibit different susceptibilities to inhibition by clavulanic acid, with KPC-2 effectively resistant to inhibition by clavulanic acid.

QM/MM simulations were carried in AMBER with the semiempirical SCC-DFTB method used for the QM region and the MM region described by AMBER’s ff14SB force field. Complete two-dimensional free energy surfaces were then obtained by QM/MM umbrella sampling (US) simulations along two reaction coordinates: one to describe the proton transfer and one to describe the nucleophilic attack involved in deacylation step of the reaction.   The weighted-histogram analysis method (WHAM) was used to generate the free energy surfaces and the minimum free energy path (MEP) on the surfaces was determined using the minimum energy pathway analysis for energy landscapes (MEPSA) software. The results show that the QM/MM screening protocol is reliable in discriminating the inhibitory activity of different covalent complexes in class A β-lactamases and can differentiate enzymes for which clavulanic acid is an effective inhibitor (TEM-1) from those for which it is not (KPC-2). Such screening techniques will prove to be useful in the development of new strategies to overcome antibiotic resistance.

Press release:

http://www.bristol.ac.uk/news/2018/june/antibiotic-arms-race.html

Full paper:

‘Multiscale simulations of clavulanate inhibition identify the reactive complex in class A β-lactamases and predict efficiency of inhibition’ by Rubén A. Fritz, Jans H. Alzate-Morales, James Spencer, Adrian J. Mulholland and Marc W. van der Kamp in Biochemistry.

Tangible physical models have an important place in the history of chemistry and biochemistry. Three-dimensional (3D) molecular models have been used as conceptual and educational tools dating back to at least von Hofmann in the 1860s. Now we are able to further our understanding by interacting with molecules in real time as they move using virtual reality (VR) technology. The multi-user VR framework has been developed by a joint team of computer science and chemistry researchers at the University of Bristol, in collaboration with developers at Interactive Scientific and Oracle Corporation, have used Oracle’s public cloud infrastructure to combine real-time molecular simulations with the latest VR technology.

The movie below shows a researcher taking hold of a fully solvated benzylpenicillin ligand and interactively guiding it to dock it within the active site of the TEM-1 β-lactamase enzyme (with both molecules fully flexible and dynamic) and generate the correct binding mode. The simulations make use of a plug-in that communicates with the GPU-accelerated molecular simulation package OpenMM via PLUMED, allowing models of the type used for conventional MD simulations to be used directly within the VR framework. Generating a benzylpenicillin docking pathway would be very difficult just using standard computer algorithms, but as can be seen in the movie, the VR framework is intuitive and easy to control, enabling the researcher to generate a physically meaningful path. Anybody wishing to try out the tasks described in the paper can download the software at https://isci.itch.io/nsb-imd, and launch their own cloud-hosted session.

Press release:

http://www.bristol.ac.uk/news/2018/july/molecular-simulations-with-vr-technology.html

Full Paper:

‘Sampling molecular conformations and dynamics in a multi-user virtual reality framework‘ by Michael O’Connor, Helen M. Deeks, Edward Dawn, Oussama Metatla, Anne Roudaut, Matthew Sutton, Becca Rose Glowacki, Rebecca Sage, Philip Tew, Mark Wonnacott, Phil Bates, Adrian J. Mulholland and David R. Glowacki in Science Advances.