The SAMPL6 SAMPLing challenge: assessing the reliability and efficiency of binding free energy calculations

Abstract

Approaches for computing small molecule binding free energies based on molecular simulations are now regularly being employed by academic and industry practitioners to study receptor-ligand systems and prioritize the synthesis of small molecules for ligand design. Given the variety of methods and implementations available, it is natural to ask how the convergence rates and final predictions of these methods compare. In this study, we describe the concept and results for the SAMPL6 SAMPLing challenge, the first challenge from the SAMPL series focusing on the assessment of convergence properties and reproducibility of binding free energy methodologies. We provided parameter files, partial charges, and multiple initial geometries for two octa-acid (OA) and one cucurbit[8]uril (CB8) host-guest systems. Participants submitted binding free energy predictions as a function of the number of force and energy evaluations for seven different alchemical and physical-pathway (i.e., potential of mean force and weighted ensemble of trajectories) methodologies implemented with the GROMACS, AMBER, NAMD, or OpenMM simulation engines. To rank the methods, we developed an efficiency statistic based on bias and variance of the free energy estimates. For the two small OA binders, the free energy estimates computed with alchemical and potential of mean force approaches show relatively similar variance and bias as a function of the number of energy/force evaluations, with the attach-pull-release (APR), GROMACS expanded ensemble, and NAMD double decoupling submissions obtaining the greatest efficiency. The differences between the methods increase when analyzing the CB8-quinine system, where both the guest size and correlation times for system dynamics are greater. For this system, nonequilibrium switching (GROMACS/NS-DS/SB) obtained the overall highest efficiency. Surprisingly, the results suggest that specifying force field parameters and partial charges is insufficient to generally ensure reproducibility, and we observe differences between seemingly converged predictions ranging approximately from 0.3 to 1.0 kcal/mol, even with almost identical simulations parameters and system setup (e.g., Lennard-Jones cutoff, ionic composition). Further work will be required to completely identify the exact source of these discrepancies. Among the conclusions emerging from the data, we found that Hamiltonian replica exchange-while displaying very small variance-can be affected by a slowly-decaying bias that depends on the initial population of the replicas, that bidirectional estimators are significantly more efficient than unidirectional estimators for nonequilibrium free energy calculations for systems considered, and that the Berendsen barostat introduces non-negligible artifacts in expanded ensemble simulations.

SAMPL7 Host-Guest Challenge Overview: Assessing the reliability of polarizable and non-polarizable methods for binding free energy calculations

“…Both analysis approaches ought to give equivalent binding free energies in certain limits, but their underlying assumptions and the amount of data available result in substantially different performance here. Here, despite its limitations in the SAMPL6 "SAMPLing" challenge [5], the Gaussian approximation was modestly superior, with, 15 of 16 binding free energies predicted within 2 kcal/mol, versus 11 of16 for FSDAM/GAFF2/OPC3_JB.…”

“…A series of density functional theory (DFT)-based methods, including ranked and non-ranked submissions, were also used here. In SAMPL6, a DFT-based approach yielded good qualitative results [4,5,43], though without geometry optimizations employed in the current challenge. Here, the QM DFT-based B2PLYPD3/SMD submissions use B3PW91 with GD3BJ [2] to treat dispersion, B2PLYPD3 for single point energy calculations [1], and the SMD implicit solvation model [90] Different submissions in this series differed in which basis set was chosen for geometry optimization.…”

“…Another important area of work is to ensure that methods which ought to be equivalent do, in fact, give equivalent results across different simulation packages. [5].…”

Automated Assessment of Binding Affinity via Alchemical Free Energy Calculations

“…Recent efforts from the field have highlighted that the results of an FEP calculation depend on several parameters other than the potential energy function, for instance input preparation, simulation conditions and choice of the free energy analysis method. 24,25 Most recently, relative binding free energy (RBFE) calculations were implemented in Cresset's structure-based drug design suite, Flare. 26,27 This approach is based on the Sire 28 and a customized version of LOMAP.…”

Capturing the Flexibility of a Protein-Ligand Complex: Binding Free Energies from Different Enhanced Sampling Techniques

“…Given this vast tool set of enhanced sampling techniques, it is a challenge of its own to select the technique that will be most efficient for a particular system of interest. For this purpose, it is helpful to know which factors affect the convergence of the respective enhanced sampling technique, which has recently been analyzed in detail for umbrella sampling 21 , and to compare the reliability and efficiency of enhanced sampling techniques, both on benchmark systems for systematic testing 22 and on demanding large-scale systems 23 . As the plasticity of Major Histocompatibility Complex I (MHC I) is the crucial factor in the intracellular selection of antigenic peptides 24 , peptide-MHC I complexes (pMHC I) are a challenging and biologically relevant test system to assess whether the enhanced sampling technique used samples the conformational dynamics sufficiently to explain biological function.…”