Alchemical Free Energy Estimators and Molecular Dynamics Engines: Accuracy, Precision and Reproducibility

The binding free energy between a ligand and its target protein is an essential quantity to know at all stages of the drug discovery pipeline. Assessing this value computationally can offer insight into where efforts should be focused in the pursuit of effective therapeutics to treat myriad diseases. In this work we examine the computation of alchemical relative binding free energies with an eye to assessing reproducibility across popular molecular dynamics packages and free energy estimators. The focus of this work is on 54 ligand transformations from a diverse set of protein targets: MCL1, PTP1B, TYK2, CDK2 and thrombin. These targets are studied with three popular molecular dynamics packages: OpenMM, NAMD2 and NAMD3. Trajectories collected with these packages are used to compare relative binding free energies calculated with thermodynamic integration and free energy perturbation methods. The resulting binding free energies show good agreement between molecular dynamics packages with an average mean unsigned error between packages of 0.5 $kcal/mol$ The correlation between packages is very good with the lowest Spearman's, Pearson's and Kendall's tau correlation coefficient between two packages being 0.91, 0.89 and 0.74 respectively. Agreement between thermodynamic integration and free energy perturbation is shown to be very good when using ensemble averaging.

Large Scale Study of Ligand-Protein Relative Binding Free Energy Calculations: Actionable Predictions from Statistically Robust Protocols

The accurate and reliable prediction of protein-ligand binding affinities can play a central role in the drug discovery process as well as in personalised medicine. Of considerable importance during lead optimisation are the alchemical free energy methods that furnish estimation of relative binding free energies (RBFE) of similar molecules. Recent advances in these methods have increased their speed, accuracy and precision. This is evident from the increasing number of retrospective as well as prospective studies employing them. However, such methods still have limited applicability in real-world scenarios due to a number of important yet unresolved issues. Here, we report the findings from a large dataset comprising over 500 ligand transformations spanning over 300 ligands binding to a diverse set of 14 different protein targets which furnish statistically robust results on the accuracy, precision and reproducibility of RBFE calculations. We use ensemble-based methods which are the only way to provide reliable uncertainty quantification given that the underlying molecular dynamics is chaotic. These are implemented using TIES (Thermodynamic Integration with Enhanced Sampling) but are equally applicable to free energy perturbation calculations for which we expect essentially very similar results. Results achieve chemical accuracy in all cases. Ensemble simulations also furnish information on the statistical distributions of the free energy calculations which exhibit non-normal behaviour. We find that the “enhanced sampling” method known as replica exchange with solute tempering degrades RBFE predictions. We also report definitively on numerous associated alchemical factors including the choice of ligand charge method, flexibility in ligand structure and the size of the alchemical region including the number of atoms involved in transforming one ligand into another. Our findings provide a key set of recommendations that should be adopted for the reliable application of RBFE methods.

The Relative Binding Position of Nck and Grb2 Adaptors Dramatically Impacts Actin-Based Motility of Vaccinia Virus

ABSTRACTPhosphotyrosine (pTyr) motifs in unstructured polypeptides orchestrate important cellular processes by engaging SH2-containing adaptors to nucleate complex signalling networks. The concept of phase separation has recently changed our appreciation of such multivalent networks, however, the role of pTyr motif positioning in their function remains to be explored. We have now explored this parameter in the assembly and operation of the signalling cascade driving actin-based motility and spread of Vaccinia virus. This network involves two pTyr motifs in the viral protein A36 that recruit the adaptors Nck and Grb2 upstream of N-WASP and Arp2/3-mediated actin polymerization. We generated synthetic networks on Vaccinia by manipulating pTyr motifs in A36 and the unrelated p14 from Orthoreovirus. In contrast to predictions, we find that only specific spatial arrangements of Grb2 and Nck binding sites result in robust N-WASP recruitment, Arp2/3 driven actin polymerization and viral spread. Our results suggest that the relative position of pTyr adaptor binding sites is optimised for signal output. This finding may explain why the relative positions of pTyr motifs are usually conserved in proteins from widely different species. It also has important implications for regulation of physiological networks, including those that undergo phase transitions.

Relative Binding Free Energy Predictions for Inhibitors of Tetrameric Influenza Virus Neuraminidase

Accurate methods to predict the free energies of protein-ligand interactions have great potential to assist rational drug design. In this work, we used molecular dynamics simulations with alchemical perturbation to predict the binding of carbohydrate-based ligands to influenza virus neuraminidase (N2). This specific drug target is a challenging test system for alchemical free energy methods because it has flexible binding site motifs. We use a molecular dynamics protocol that works for longer time scales than are often reported in previous molecular dynamics studies of N2. We demonstrated that N2-ligand complex stability and that accurate N2 150-loop dynamics, on a 350 ns time scale, are both force field-dependent (AMBER99SB-ILDN, GAFF and TIP4P water are required). Further, we showed that crystallographic waters must be retained to maintain stability of N2-ligand complexes over 350 ns. Using our modelling protocol, we were able to determine relative binding free energy values between neuraminidase ligands which correlated strongly with experimental differences in pIC50 values (R = -0.96, ρ = 0.86, N = 13, sig < 0.0001). It is anticipated that the molecular dynamics protocol and the relative binding free energy methods reported here, will both be useful in expediting the discovery of novel therapeutics for N2 and other homologous glycohydrolases.