Free-Energy Workshop Abstract Page

Session 1: Presentation 1

Chris Oostenbrink, BOKU, Vienna, Austria

Robust free-energy calculations from a small number of simulations

In this presentation I will demonstrate our efforts to compute alchemical free-energy differences following robust statistical mechanics, but using a limited number of simulations. In a first approach [1] , we combine approximate the nonpolar contribution to the free energy from a one-step perturbation from a single reference state, while the electrostatic contributions are obtained from simulations at the endstates and a third-power polynomial fit to the derivatives required in thermodynamic integration. We apply this approach to all amino-acid sidechains to compute thefree energy of solvation of the sidechain analogs and the freeenergy of mutation in small model peptides. In a second approach [2], we use a reweighting scheme to correctly interpolate an entire thermodynamic integration curve, involving any kind of interaction, from a small set of simulations. Using this extended TI, one immediately gets more precise estimates of the free energydifferences and can very efficiently decide where to focus additional simulation efforts.

1. A. de Ruiter and C. Oostenbrink J. Chem. Theory Comp. 8 (2012) 3686 – 3695
2. A. de Ruiter and C. Oostenbrink J. Chem. Theory Comput. 12 (2016) 4476 – 4486

Session 1: Presentation 2

Berk Hess, KTH, Stockholm, Sweden

Optimal sampling along reaction coordinates

Sampling of conformational transitions in (bio-)molecules decreases exponentially with the high of barriers in the free-energy landscape. Therefore exponential acceleration of sampling can be achieved by applying methods that flatten the sampled distribution. There are several methods available that do this, for instance metadynamics and the accelerated weight histogram method, which will be presented here. But when modifying sampling from the original Boltzmann sampling, the issue of a metric arises. One can make the sampling along a reaction coordinate uniform, but then a non-linear transformation of that reaction coordinate will result in non uniform sampling. Here I present a metric that measures the friction along one or multidimensional reaction coordinates. This provides useful insight into difficult regions along transitions and enables optimization of the sampling.

Session 1: Presentation 3

Peter Coveney, UCL

Ensemble-based molecular dynamics: principles and applications

In MD simulations, macroscopic properties corresponding to experimental observables are defined in terms of ensemble averages. Foundational aspects of statistical mechanics are, however, largely overlooked in most textbooks and research articles that purport to compute macroscopic behaviour from microscopic descriptions based on classical mechanics. Ensemble-averaged behaviour needs to be computed in a manner which is truly consistent with the methods of statistical mechanics. In the past few years, we have constructed ensemble-based molecular dynamics approaches which provide a route to reliable predictions of free energies that meet the requirements of speed, accuracy, precision, and reproducibility. Our predictions using these approaches, some of them performed blindly, are in good agreement with experimental findings.

Session 1: Presentation 4

Davide Branduardi, Schrödinger, Inc., USA

Schrödinger FEP+: calculating relative binding free energiesin lead optimization

The prediction of relative free energies of binding within a congeneric ligand series has been considered a potential key contribution to the process of ligand optimization since the advent of computer aided drug discovery. Whenever accurately computed, free energy predictions would allow reducing the effort spent in compound synthesis and decrease compound design cycle times. Traditionally, fast empirical approaches like rigid-receptor docking dominated the field for decades due to their low computational cost but have proven to be very limited in their predictive power in lead optimization and often require a large amount of experience to interpret as results can be highly target specific. Nowadays we are witnessing the beginning of a new era for relative binding free energy predictions due to the emergence of freeenergy perturbation (FEP) calculations. Although this technique was established more than twenty years ago, it is now seeing a broader application thanks to a number of technological advances.

In this talk I’ll introduce the main features of Schrodinger’s FEP+ solution, a ligand optimization protocol which combines the recent advances in the field. FEP+ features high performance MD code that makes use of GPU computing, augmented by enhanced sampling schemes such as “replica-exchange solute tempering”, so to achieve sufficient effective simulation time to provide converged free energyestimates with modest hardware investment. An error analysis based on cycle-closure correction provides a measure of reliability of the calculations. Additionally FEP+ includes the accurate all-atom forcefield OPLS3 and automated setup-protocols which enable free energycalculations to be performed in an industry environment. Beyond the technical overview, I will present a few cases from literature and internal efforts where FEP+ has been used retrospectively and prospectively in order to outline Schrödinger’s efforts in extending the domain of applicability of this technique.

Session 2: Presentation 1

Antonia Mey, University of Edinburgh, United Kingdom

Current limits of binding free energy calculations

The basic principles of alchemical free energy calculations are well-understood, but reduction to practice still requires uncertain decisions that critically affect final outcome. This talk will summarise lessons learned from applications of an alchemical free energy calculation workflow developed in my group over the last 5 years.

The underlying science mandates that excess free energy changes are uniquely encoded by a potential energy function, given sufficient configurational averaging. In practice code-specific implementation details may inadvertently prevent internal consistency of computed freeenergies, let alone reproducibility with other software packages.

In a retrospective application to thrombin ligands, remarkable gains in accuracy over docking methods are observed, but it is difficult to distil numbers into simple ligand design guidelines. For RSL lectin binding ligands the generic AM1-BCC/GAFF force-field is on par with the best docking methods we tested, but the specialized GLYCAM force-field significantly improves results.

In blind prediction studies significant correlations with experiment were achieved for three classes of host-guest systems as part of the last SAMPL5 competitions. However systematic deviations with experimental data indicate that the underlying physics is incompletely captured by the models. Blind predictions of affinities for HSP90 and FXR ligands as part of D3R grand challenges did not achieve high discriminatory power, despite encouraging retrospective results on test sets. Post-mortem analyses suggest that uncertainties in binding site water placement and ligand binding modes, inaccurate potential energyfunctions, and approximate treatment of electrostatic interactions, can critically cloud the validity of computed binding free energies.

Session 2: Presentation 2

Hannah Bruce Macdonald, University of Southampton

Predicting water networks and relative ligand binding freeenergies in proteins using grand canonical Monte Carlo

Understanding the location and energetics of water molecules is needed for rational drug design. We present our grand canonical Monte Carlo (GCMC) method for predicting water locations and binding freeenergies in proteins.[1] Further, we can calculate binding free energies of ligands with different active-site water networks by combining thermodynamic integration and GCMC in a single simulation (GCTI).

The energetics of an active site water molecule can be vital to drug design, through either stabilising the bound complex, or releasing entropy upon displacement. Rationalising how water molecules should be treated is difficult – whether a water molecule should be retained in a bound structure, or if it should be displaced to recover entropy and allow for direct protein-ligand interactions, is unclear.[2]

Methods exist to calculate the binding free energies of water molecules through decoupling from the surrounding environment. This requires a priori knowledge of the water location, and can only calculate the freeenergy of a single water molecule per simulation. GCMC is able to calculate the free energy of multiple water networks in a single simulation, and as the water location is predicted automatically as part of the simulation, no experimental knowledge of hydration site location is required. Free energies calculated using GCMC are consistent with the results of decoupling simulations.

While it is important to understand the hydration thermodynamics in protein systems, the most important property remains the ligand binding free energy. While conventional free energy simulations can calculate relative binding free energies of ligands, issues can arise when the water network changes with the ligand, particularly in the case of occluded cavities.[3] We show that GCMC can now be performed in conjunction with ligand binding free energy simulations. This allows for direct calculation of relative free energies between ligands with different water occupancies within a single simulation. We refer to this combination of GCMC and thermodynamic integration as GCTI. Experimental data for Scytalone Dehydratase for congenic ligands, where displacement of a water molecule results in a large increase in affinity, have been reproduced using the GCTI method. This powerful GCTI method is therefore able to directly provide the free energy difference between ligands while automatically accounting for hydration effects.

1. G.A. Ross, M.S. Bodnarchuk, J.W. Essex, J. Am. Chem. Soc. 173 (2015) 14930-14943.
2. C. Barillari, J. Taylor, R. Viner, J. W. Essex, J. Am. Chem. Soc. 129 (2007) 2577-2587.
3. J. Michel, J. Tirado-Rives, W. L. Jorgensen, J. Am. Chem. Soc. 131 (2009) 15403-15411.

Session 2: Presentation 3

Phil Biggin, Oxford

Predictions of ligand selectivity from absolute binding free energy calculations

Binding selectivity is an essential property for effective chemical probes used in preclinical target validation and is a requirement for the development of a safe drug. The rational design of selectivity adds considerable complexity to an already complex problem. Computationally, the prediction of binding selectivity is a challenge and general applicable methodologies are still not available to the medicinal chemistry community. Absolute binding free energy calculations based on alchemical pathways provide a rigorous framework for affinity predictions and could thus offer a general approach to the problem. As there has been much recent progress in the methodological aspects with respect to this approach, we have interested in assessing its performance when utilizing all of advice falling under the umbrella of “best practice”. We have evaluated the performance of free energy calculations based on molecular dynamics for the prediction of selectivity by estimating the affinity profile of three bromodomain inhibitors across multiple bromodomain families, and by comparing the results to isothermal titration calorimetry data. We have considered two case studies: In the first case, the affinities of a broad-spectrum inhibitor for 22 bromodomains were calculated, and returned an accuracy (mean unsigned error of 1.76 kcal/mol and Pearson correlation of 0.48). In the second case, the affinities of two similar ligands for seven bromodomains were calculated and returned excellent agreement with experiment (mean unsigned error of 0.81 kcal/mol and Pearson correlation of 0.75). In this test case, we also show how the preferred binding orientation of a ligand for different proteins can be determined via free energy calculations – a property which is inherently problematic for relative free energy calculation methods. In this talk I will outline these results some of which were recently published.[1]

1. Aldeghi M, Heifetz A, Bodkin MJ, Knapp S, Biggin PC (2017) Predictions of ligand selectivity from absolute binding free energycalculations. J Am Chem Soc 139:946-957. doi:10.1021/jacs.6b11467

Session 2: Presentation 4

Alexander Heifetz, Evotec

Accurate Assessment of Protein-Ligand Interaction Energy in Seconds with Quantum Mechanics

The understanding of binding interactions between any protein and a small molecule is a cornerstone of any efficient structure-based drug design (SBDD) process. X-ray crystallography and homology modelling are the main source of structural information required for rational SBDD. However, even with the crystal structure in hand, “visual inspection” and force field-based molecular mechanics calculations often used for the rationalization of ligand-protein potency cannot always explain the full complexity of the molecular interactions. Quantum mechanical (QM) approach was always considered as promising direction to achieve this goal however, traditional QM are not feasible for large biological systems, due to their high computational cost.

FMO method offers a considerable computational speed-up over traditional QM methods. One of the key features of the FMO approach is that it can provide a list of the interactions formed between the ligand and the receptor and a chemically intuitive breakdown of these interactions. Such information is essential for medicinal chemists to be able to rationally approach modification of lead compounds in order to increase favourable interactions. We will demonstrate the prospective application of FMO method in drug-discovery programs.

Recently, we have demonstrated that FMO can be even faster (secs instead of hours) without compromising the accuracy by combining it with density-functional tight-binding (DFTB) method. We will exemplified the use of FMO-DFTB in three GPCR-ligand systems.

Session 3: Presentation 1

Hannes H Loeffler, SCD, STFC Daresbury, UK

FESetup: automating relative alchemical free energy simulations

The alchemical free energy (AFE) method has been applied in various forms since the early days of using computers for molecular simulation. Currently there is renewed interest in this simulation method as advances in computing technology and simulation theory make it increasingly possible to compute free energies for large systems of biological interest. As AFE is dependent on sufficient sampling it is, computationally, a relatively expensive way of computing free energies.

AFE is rooted in a rigorous treatment of statistical thermodynamics and as such is principally capable of producing accurate results. Sampling and quality of force files are fundamental issues but a more basic problem is that the method is often difficult to use with present day simulation packages. The setup can be very time–consuming, laborious and very error–prone. Compared to standard MD simulation, AFE calculation is more difficult to carry out than it should be.

I will present a tool, FESetup[1], which can help in setting up relative alchemical free energy simulations. FESetup is developed as part of the software support effort within the Collaborative Computational Project for Biomolecular Simulation (CCPBioSim) [2]. FESetup automates the writing of perturbed topology files as much as possible by using a maximum common substructure search (MCSS) strategy to exploit topological similarity in a pair of ligand molecules. The aim of the software is to minimise the human bottleneck by helping the researcher to focus more on scientific problems and much less on the intricacies of a particular software. Based on the challenges that I met in developing this tool, I will discuss what is required to ensure a consistent simulation setup, where improvements are still needed and what the future opportunities for the field are. Results from an ongoing reproducibility study will be presented too.

1. DOI:10.1021/acs.jcim.5b00368, http://www.ccpbiosim.ac.uk/fesetup/
2. http://www.ccpbiosim.ac.uk/

Session 3: Presentation 2

Gianni De Fabritiis, UPF, Barcelona, Spain

HTMD: High-throughput molecular dynamics for molecular discovery

Recent advances in molecular simulations have allowed scientists to investigate slower biological processes than ever before. Together with these advances came an explosion of data that has transformed a
traditionally computing-bound into a data-bound problem. Here, we present HTMD, a programmable, extensible platform written in Python that aims to solve the data generation and analysis problem as well as increase reproducibility by providing a complete workspace for simulation-based discovery. So far, HTMD includes system building for CHARMM and AMBER force fields, projection methods, clustering, molecular simulation production, adaptive sampling, an Amazon cloud interface, Markov state models, and visualization. As a result, a single, short HTMD script can lead from a PDB structure to useful
quantities such as relaxation time scales, equilibrium populations,metastable conformations, and kinetic rates.

Session 3: Presentation 3

Erik Lindahl, KTH, Stockholm, Sweden

Efficient high-throughput free energy calculation and parameterization with GROMACS & STaGE

GROMACS provides a large number of options for free energycalculations, but they are not all created equal. Here, I will present some of the recent possibilities that enable you to perform both more accurate and more efficient free energy calculations than what has been possible before, and how this can be combined with important features such as acceleration using CUDA, Intel Xeon Phi, or OpenCL-compatible cards. I will also present some of the tools we have developed to perform fully automatic parameterization of arbitrary small ligands for all major force fields, and distributed computing techniques that both automate large-scale free energy screening and provide new options to optimize what calculations to perform.

Session 3: Presentation 4

Bert de Groot, MPI, Goettingen, Germany

Large-scale stability and affinity estimates through alchemical freeenergy calculations using pmx.

The pmx framework allows for the generation of fully automated, high quality hybrid GROMACS topologies for alchemical freeenergy calculations. I will present a large-scale amino acid mutation scan of protein stability using several contemporary molecular mechanics force fields. Comparison to ITC data provides a thermodynamic benchmark of these force fields and combining the results from multiple force fields in a consensus approach is shown to outperform the accuracy of individual force fields. In addition, support for nucleic acid substitutions and ligand modification has been added, for which first results will be presented.

Session 4: Presentation 1

Philip Fowler, John Ratcliffe Hospital, University of Oxford

De novo prediction of resistance to trimethoprim in Staphylococcus aureus.

The discovery of antibiotics was one of humanity’s greatest achievements in the twentieth century, however, the evolution of antibiotic resistance by pathogens now threatens many advances of modern medicine. There is an urgent need for improved diagnostic tools so that resistant infections can be identified and treated appropriately. Analysis of whole-genome sequence data generated on affordable high-throughput platforms has the potential to allow resistant infections to be more rapidly and cheaply diagnosed in the clinic than conventional culture based approaches. A key limitation of this approach is that it cannot identify rare or previously unseen mutations. I shall demonstrate that a well-established class of methods from computational chemistry, alchemical free energy methods, can successfully predict the resistance phenotype of a series of mutations in Staphylococcus aureus dihydrofolate reductase identified by whole-genome sequencing of patient infections. Not only can the method predict which mutations cause resistance (and which do not), but it also can predict the minimum inhibitory concentration of trimethoprim for each mutation with reasonable accuracy. Having established that the approach has the potential to be successful, I shall briefly discuss how this and similar approaches are being used in a major international tuberculosis collaboration, CRyPTIC, that was launched in March 2016. The goal of this project is to collect
100,000 M. tuberculosis samples from across Africa, Asia, Europe and the Americas. Each sample will have both its drug susceptibility profile established and its whole genome sequenced.

Session 4: Presentation 2

Daniel Seeliger, Boehringer Ingelheim Pharma Biberach, Ingelheim am Rhein, Germany

Towards free-energy calculations as routine applications in industrial drug discovery

Opportunities for the development of new therapeutics has expanded from small-molecule approaches towards novel therapeutic modalities ranging from peptides, protein-based therapeutics, RNA, cancer vaccines, CAR-T cells, and gene therapy. While research budgets are generally not substantially increasing, opportunities for investments in therapeutics research have diversified considerably. For classical small-molecule drug discovery these developments denote increasing pressure to become more efficient and to focus on more sophisticated design and prediction of molecule properties. The presentation will outline the current state of computational design of small molecules at Boehringer Ingelheim and discuss current gaps and shortcomings.

Session 4: Presentation 3

Miguel Machuqueiro, University of Lisbon

pKa shifts in protein and membrane binding

Pedro Reis, Tomás Silva, Bruno Victor, Diogo Vila Viçosa and Miguel Machuqueiro
Centro de Química e Bioquímica, Faculdade de Ciências, Universidade de Lisboa, 1749-016 Lisboa, Portugal; Email: machuque@ciencias.ulisboa.pt
pH is a crucial physicochemical property that affects most biomolecules. Changes in protonation equilibrium of susceptible sites will modify the electrostatic environment and, consequently, have an effect on the molecular structure, stability and catalysis.[1] The pKa values of common titrable sites in peptides or other simple organic molecules can be significantly influenced by changes in solvent mixture, by direct interaction with a protein bindig pocket or due to insertion in a lipid bilayer.[2-3] In this work, we present our latest results on the pKa calculations of peptides at the water/membrane interface and known acetylcholinesterase (AChE) inhibitors bound to their receptor. We take advantage of the recent extensions to the CpHMD-L methodology [4] and apply it to these different systems, namely, the model Ala-based pentapeptides that have already been well characterized in water by Pace and co-workers,[5] and to donepezil and galantamine, two commercially available drugs that are inhibitors of AChE.[6] We acknowledge financial support from FCT through projects PTDC/QEQCOM/5904/2014, UID/MULTI/00612/2013.

1. J. Matthew, Annu. Rev. Biophys. Biophys. Chem. 14 (1985) 387.
2. J. F. Hunt, P. Rath, K. J. Rothschild, D. M. Engelman, Biochemistry 36 (1997) 15177.
3. P. R. Magalhães, M. Machuqueiro, A. M. Baptista, Biophys. J. 108 (2015) 2282.
4. D. Vila-Viçosa, V. H. Teixeira, A. M. Baptista, M. Machuqueiro, JCTC 11 (2015) 2367.
5. R. L. Thurlkill, G. R. Grimsley, J. M. Scholtz, C. N. Pace, Protein Sci. 15 (2006) 1214.
6. P. Draczkowski, A. Tomaszuk, P. Halczuk, M. Strzemski, D. Matosiuk, K. Jozwiak, Biochim. Biophys. Acta, 1860 (2016) 967.

Session 4: Presentation 4

Gary Tresadern, Janssen, Antwerp, Belgium

Free Energy Perturbation Applied in Drug Discovery at Janssen: Learnings so far…

Accurate prediction of binding energies between small molecules and their protein targets has the potential to make a massive impact on drug discovery. Free energy perturbation is a well-established theoretical approach to calculate the relative binding energy difference between two molecules. By perturbing one molecule to the other in solvent and in protein, then taking the difference, the computationally demanding step of calculating absolute binding energies is circumvented. In recent years, free energy perturbation has benefited from various practical and methodological improvements that suggest its suitability for large scale implementation in industrial drug discovery.[1],[2] Therefore, we have been assessing the performance of freeenergy perturbation in various projects at Janssen.[3],[4],[5] In this talk we will review some of our applications, discuss the successes and failures, and summarize our current status.

1. Wang, L.; et al. Accurate and Reliable Prediction of Relative Ligand Binding Potency in Prospective Drug Discovery by Way of a Modern Free-Energy Calculation Protocol and Force Field. J. Am.Chem. Soc. 2015, 137, 2695−2703.
2. Harder, E.; et al OPLS3: A Force Field Providing Broad Coverage of Drug-like Small Molecules and Proteins. J. Chem. Theory Comput. 2016, 12, 281−296.
3. Rombouts, F. J. R.; et al. Pyrido[4,3-e][1,2,4]triazolo[4,3-a]pyrazines as Selective, Brain Penetrant Phosphodiesterase 2 (PDE2) Inhibitors. ACS Med. Chem. Lett. 2015, 6, 282−286.
4. Ciordia, M.; et al Application of Free Energy Perturbation for the Design of BACE1 Inhibitors. J. Chem. Inf. Model. 2016, 56, 1856−1871.
5. Keränen, H.; et al Acylguanidine Beta Secretase 1 Inhibitors: A Combined Experimental and Free Energy Perturbation Study. J. Chem. Theory Comput. 2017, Article ASAP.