Standard Binding Free-Energy Calculations: How Far Are We from Automation?

Identification of novel inhibitors for epidermal growth factor receptor tyrosine kinase using absolute binding free-energy simulations

Mutations in the kinase domain of the epidermal growth factor receptor (EGFR), a critical biological macromolecule involved in cell growth and division, can lead to drug resistance in patients undergoing chemotherapy with kinase inhibitors. Notably, the emergence of the C797S mutation poses new challenges for targeted EGFR therapy, highlighting the urgent need for agents effective against this triple mutation (L858R/T790M/C797S, EGFR™). Building on our previous finding that sulfonyl and piperidinyl groups significantly contribute to the EGFR™-inhibitor interactions, we have identified the best-in-class inhibitors containing these groups through functional-group-based screening and formally exact absolute binding free-energy calculations. Our new strategy offers greater flexibility than traditional workflows leaning on relative binding free-energy calculations and accommodates ligands with substantial structural variations. The result shows that the top candidate exhibits a binding affinity of -15.8 kcal/mol towards the EGFR™ mutant, surpassing BLU-945, a state-of-the-art fourth-generation inhibitor with a binding free energy of -12.6 kcal/mol. Subsequent free-energy decomposition indicates that the presented top candidate primarily enhances interactions with the K745, D800 and R841 residues, suggesting its potential to overcome resistance from the C797S mutation. Notably, K745 forms highly favorable hydrogen bonds and cation-π interactions with C6. Targeting lysine has emerged as a promising strategy, especially in cases where the C797S mutation renders traditional covalent inhibitors ineffective. We propose that these novel inhibitors represent promising drug candidates for non-small cell lung cancer treatment and offer new strategies to overcome drug resistance caused by EGFR mutation.

DeepConf: Leveraging ANI-ML Potentials for Exploring Local Minima with Application to Bioactive Conformations

Here, we introduce a low energy conformer generation algorithm using ANI-ML potentials at the DFT accuracy and benchmark in reproducing bioactive conformations. We show that the method is efficient when the initial starting structure is far from equilibrium, when the ML potentials are stuck in nonsmooth regions, and when the quality of the conformers in a less conformer size is demanded. We specifically focus on conformations due to rotations around the single bonds. For the first time, we assess the performance of ANI-ML potentials using our conformer generation algorithm, DeepConf, in addition to previously reported Auto3D (J. Chem. Inf. Model. 2022, 62, 5373-5382) using the same potentials to reproduce bioactive conformations as well as providing a guideline for bioactive conformation evaluation processes. Our results show that the ANI-ML potentials can reproduce the bioactive conformations with mean value of the root-mean-square-deviation (RMSD) less than 0.5 Å, outperforming the limit of conventional methods. The code offers several features including but not limited to geometry optimization, fast conformer generations via single point energies (SPE), different minimization algorithms, different ML-potentials, or high-quality conformers in the smallest amount of ensemble sizes. It is available free of charge (documentation and test files) at https://github.com/otayfuroglu/DeepConf.

Large-Scale AI-Based Structure and Activity Prediction Analysis of ShK Domain Peptides from Sea Anemones in the South China Sea

Sea anemone peptides represent a valuable class of biomolecules in the marine toxin library due to their various structures and functions. Among these, ShK domain peptides are particularly notable for their selective inhibition of the Kv1.3 channel, holding great potential for applications in immune regulation and the treatment of metabolic disorders. However, these peptides' structural complexity and diversity have posed challenges for functional prediction. In this study, we compared 36 ShK domain peptides from four species of sea anemone in the South China Sea and explored their binding ability with Kv1.3 channels by combining molecular docking and dynamics simulation studies. Our findings highlight that variations in loop length, residue composition, and charge distribution among ShK domain peptides affect their binding stability and specificity. This work presents an efficient strategy for large-scale peptide structure prediction and activity screening, providing a valuable foundation for future pharmacological research.

Flexible framework of computing binding free energy using the energy representation theory of solution

Host-guest binding plays a crucial role in the functionality of various systems, and its efficiency is often quantified using the binding free energy, which represents the free-energy difference between the bound and dissociated states. Here, we propose a methodology to compute the binding free energy based on the energy representation (ER) theory of solution, which enables us to evaluate the free-energy difference between the systems of interest with the molecular dynamics (MD) simulations. Unlike the other free-energy methods, such as the Bennett acceptance ratio (BAR), the ER theory does not require the MD simulations for hypothetical intermediate states connecting the systems of interest, leading to reduced computational costs. By constructing the thermodynamic cycle of the binding process that is suitable for the ER theory, a robust calculation of the binding free energy is realized. We apply the present method to the self-association of N-methylacetamide in different solvents and the binding of aspirin to β-cyclodextrin (CD) in water. In the former case, the present method estimates that the binding free energy decreases as the solvent polarity decreases. This trend is consistent with the experimental finding. For the latter system, the binding free energies for the two representative CD-aspirin bound complexes, primary (P) and secondary (S) complexes, are estimated to be -5.2 ± 0.1 and -5.03 ± 0.09 kcal mol-1, respectively. These values are satisfactorily close to those from the BAR method [-4.2 ± 0.2 and -4.1 ± 0.2 kcal mol-1 for P and S, respectively]. Furthermore, the interaction-energy component analysis reveals that the van der Waals interaction between aspirin and CD dominantly contributes to the stabilization of the bound complexes, which is in harmony with the well-known binding mechanism in the CD systems.

Overcoming Sampling Issues and Improving Computational Efficiency in Collective-Variable-Based Enhanced-Sampling Simulations: A Tutorial

This tutorial is designed to help users overcome sampling challenges and improve computational efficiency in collective-variable (CV)-based enhanced-sampling, or importance-sampling, simulations. Toward this end, we introduce well-tempered metadynamics-extended adaptive biasing force (WTM-eABF) and its integration with Gaussian accelerated molecular dynamics (GaMD). Additionally, use will be made of a method for identifying the least-free-energy pathway (LFEP) and multiple concurrent pathways on high-dimensional free-energy surfaces. We illustrate these sampling techniques with the conformational equilibria of trialanine and chignolin in aqueous solution as test cases. This tutorial assumes that the user has prior experience with molecular dynamics (MD) simulations, in general, with the popular program NAMD, and to some extent with Colvars, the module for CV-based calculations. This tutorial can, however, in large measure be used in conjunction with alternate MD engines that support the Colvars module such as GROMACS, LAMMPS, and Tinker-HP.

Automated Adaptive Absolute Binding Free Energy Calculations

Alchemical absolute binding free energy (ABFE) calculations have substantial potential in drug discovery, but are often prohibitively computationally expensive. To unlock their potential, efficient automated ABFE workflows are required to reduce both computational cost and human intervention. We present a fully automated ABFE workflow based on the automated selection of λ windows, the ensemble-based detection of equilibration, and the adaptive allocation of sampling time based on inter-replicate statistics. We find that the automated selection of intermediate states with consistent overlap is rapid, robust, and simple to implement. Robust detection of equilibration is achieved with a paired t-test between the free energy estimates at initial and final portions of a an ensemble of runs. We determine reasonable default parameters for all algorithms and show that the full workflow produces equivalent results to a nonadaptive scheme over a variety of test systems, while often accelerating equilibration. Our complete workflow is implemented in the open-source package A3FE (https://github.com/michellab/a3fe).

Do All Paths Lead to Rome? How Reliable is Umbrella Sampling Along a Single Path?

Molecular dynamics (MD) simulations are widely applied to estimate absolute binding free energies of protein-ligand and protein-protein complexes. A routinely used method for binding free energy calculations with MD is umbrella sampling (US), which calculates the potential of mean force (PMF) along a single reaction coordinate. Surprisingly, in spite of its widespread use, few validation studies have focused on the convergence of the free energy computed along a single path for specific cases, not addressing the reproducibility of such calculations in general. In this work, we therefore investigate the reproducibility and convergence of US along a standard distance-based reaction coordinate for various protein-protein and protein-ligand complexes, following commonly used guidelines for the setup. We show that repeating the complete US workflow can lead to differences of 2-20 kcal/mol in computed binding free energies. We attribute those discrepancies to small differences in the binding pathways. While these differences are unavoidable in the established US protocol, the popularity of the latter could hint at a lack of awareness of such reproducibility problems. To test if the convergence of PMF profiles can be improved if multiple pathways are sampled simultaneously, we performed additional simulations with an adaptive-biasing method, here the accelerated weight histogram (AWH) approach. Indeed, the PMFs obtained from AHW simulations are consistent and reproducible for the systems tested. To the best of our knowledge, our work is the first to attempt a systematic assessment of the pitfalls in one the most widely used protocols for computing binding affinities. We anticipate therefore that our results will provide an incentive for a critical reassessment of the validity of PMFs computed with US, and make a strong case to further benchmark the performance of adaptive-biasing methods for computing binding affinities.

Understanding the impact of binding free energy and kinetics calculations in modern drug discovery

INTRODUCTION

For rational drug design, it is crucial to understand the receptor-drug binding processes and mechanisms. A new era for the use of computer simulations in predicting drug-receptor interactions at an atomic level has begun with remarkable advances in supercomputing and methodological breakthroughs.

AREAS COVERED

End-point free energy calculation methods such as Molecular Mechanics/Poisson Boltzmann Surface Area (MM/PBSA) or Molecular-Mechanics/Generalized Born Surface Area (MM/GBSA), free energy perturbation (FEP), and thermodynamic integration (TI) are commonly used for binding free energy calculations in drug discovery. In addition, kinetic dissociation and association rate constants (k off and k on ) play critical roles in the function of drugs. Nowadays, Molecular Dynamics (MD) and enhanced sampling simulations are increasingly being used in drug discovery. Here, the authors provide a review of the computational techniques used in drug binding free energy and kinetics calculations.

EXPERT OPINION

The applications of computational methods in drug discovery and design are expanding, thanks to improved predictions of the binding free energy and kinetic rates of drug molecules. Recent microsecond-timescale enhanced sampling simulations have made it possible to accurately capture repetitive ligand binding and dissociation, facilitating more efficient and accurate calculations of ligand binding free energy and kinetics.

Understanding the Reversible Binding of a Multichain Protein-Protein Complex through Free-Energy Calculations

We demonstrate that the binding affinity of a multichain protein-protein complex, insulin dimer, can be accurately predicted using a streamlined route of standard binding free-energy calculations. We find that chains A and C, which do not interact directly during binding, stabilize the insulin monomer structures and reduce the binding affinity of the two monomers, therefore enabling their reversible association. Notably, we confirm that although classical methods can estimate the binding affinity of the insulin dimer, conventional molecular dynamics, enhanced sampling algorithms, and classical geometrical routes of binding free-energy calculations may not fully capture certain aspects of the role played by the noninteracting chains in the binding dynamics. Therefore, this study not only elucidates the role of noninteracting chains in the reversible binding of the insulin dimer but also offers a methodological guide for investigating the reversible binding of multichain protein-protein complexes utilizing streamlined free-energy calculations.

Binding Thermodynamics of Fourth-Generation EGFR Inhibitors Revealed by Absolute Binding Free Energy Calculations

The overexpression or mutation of the kinase domain of the epidermal growth factor receptor (EGFR) is strongly associated with non-small-cell lung cancer (NSCLC). EGFR tyrosine kinase inhibitors (TKIs) have proven to be effective in treating NSCLC patients. However, EGFR mutations can result in drug resistance. To elucidate the mechanisms underlying this resistance and inform future drug development, we examined the binding affinities of BLU-945, a recently reported fourth-generation TKI, to wild-type EGFR (EGFRWT) and its double-mutant (L858R/T790M; EGFRDM) and triple-mutant (L858R/T790M/C797S; EGFRTM) forms. We compared the binding affinities of BLU-945, BLU-945 analogues, CH7233163 (another fourth-generation TKI), and erlotinib (a first-generation TKI) using absolute binding free energy calculations. Our findings reveal that BLU-945 and CH7233163 exhibit binding affinities to both EGFRDM and EGFRTM stronger than those of erlotinib, corroborating experimental data. We identified K745 and T854 as the key residues in the binding of fourth-generation EGFR TKIs. Electrostatic forces were the predominant driving force for the binding of fourth-generation TKIs to EGFR mutants. Furthermore, we discovered that the incorporation of piperidinol and sulfone groups in BLU-945 substantially enhanced its binding capacity to EGFR mutants. Our study offers valuable theoretical insights for optimizing fourth-generation EGFR TKIs.