Accurate PDZ/Peptide Binding Specificity with Additive and Polarizable Free Energy Simulations

System-Specific Parameter Optimization for Nonpolarizable and Polarizable Force Fields

The accuracy of classical force fields (FFs) has been shown to be limited for the simulation of cation-protein systems despite their importance in understanding the processes of life. Improvements can result from optimizing the parameters of classical FFs or by extending the FF formulation by terms describing charge transfer (CT) and polarization (POL) effects. In this work, we introduce our implementation of the CTPOL model in OpenMM, which extends the classical additive FF formula by adding CT and POL. Furthermore, we present an open-source parametrization tool, called FFAFFURR, that enables the (system-specific) parametrization of OPLS-AA and CTPOL models. The performance of our workflow was evaluated by its ability to reproduce quantum chemistry energies and by molecular dynamics simulations of a zinc-finger protein.

What to Make of Zero: Resolving the Statistical Noise from Conformational Reorganization in Alchemical Binding Free Energy Estimates with Metadynamics Sampling

We introduce the self-relative binding free energy (self-RBFE) approach to evaluate the intrinsic statistical variance of dual-topology alchemical binding free energy estimators. The self-RBFE is the relative binding free energy between a ligand and a copy of the same ligand, and its true value is zero. Nevertheless, because the two copies of the ligand move independently, the self-RBFE value produced by a finite-length simulation fluctuates and can be used to measure the variance of the model. The results of this validation provide evidence that a significant fraction of the errors observed in benchmark studies reflect the statistical fluctuations of unconverged estimates rather than the models' accuracy. Furthermore, we find that ligand reorganization is a significant contributing factor to the statistical variance of binding free energy estimates and that metadynamics-accelerated conformational sampling of the torsional degrees of freedom of the ligand can drastically reduce the time to convergence.

Performance and Analysis of the Alchemical Transfer Method for Binding-Free-Energy Predictions of Diverse Ligands

The Alchemical Transfer Method (ATM) is herein validated against the relative binding-free energies (RBFEs) of a diverse set of protein-ligand complexes. We employed a streamlined setup workflow, a bespoke force field, and AToM-OpenMM software to compute the RBFEs of the benchmark set prepared by Schindler and collaborators at Merck KGaA. This benchmark set includes examples of standard small R-group ligand modifications as well as more challenging scenarios, such as large R-group changes, scaffold hopping, formal charge changes, and charge-shifting transformations. The novel coordinate perturbation scheme and a dual-topology approach of ATM address some of the challenges of single-topology alchemical RBFE methods. Specifically, ATM eliminates the need for splitting electrostatic and Lennard-Jones interactions, atom mapping, defining ligand regions, and postcorrections for charge-changing perturbations. Thus, ATM is simpler and more broadly applicable than conventional alchemical methods, especially for scaffold-hopping and charge-changing transformations. Here, we performed well over 500 RBFE calculations for eight protein targets and found that ATM achieves accuracy comparable to that of existing state-of-the-art methods, albeit with larger statistical fluctuations. We discuss insights into the specific strengths and weaknesses of the ATM method that will inform future deployments. This study confirms that ATM can be applied as a production tool for RBFE predictions across a wide range of perturbation types within a unified, open-source framework.

Structure and Dynamics of Confined Water Inside Diphenylalanine Peptide Nanotubes

Diphenylalanine (FF) peptides exhibit a unique ability to self-assemble into nanotubes with confined water molecules playing pivotal roles in their structure and function. This study investigates the structure and dynamics of diphenylalanine peptide nanotubes (FFPNTs) using all-atom molecular dynamics (MD) and grand canonical Monte Carlo combined with MD (GCMC/MD) simulations with both the CHARMM additive and Drude polarizable force fields. The occupancy and dynamics of confined water molecules were also examined. It was found that less than 2 confined water molecules per FF help stabilize the FFPNTs on the x-y plane. Analyses of the kinetics of confined water molecules revealed distinctive transport behaviors for bound and free water, and their respective diffusion coefficients were compared. Our results validate the importance of polarizable force field models in studying peptide nanotubes and provide insights into our understanding of nanoconfined water.

EnCPdock: a web-interface for direct conjoint comparative analyses of complementarity and binding energetics in inter-protein associations

CONTEXT

Protein-protein interaction (PPI) is a key component linked to virtually all cellular processes. Be it an enzyme catalysis ('classic type functions' of proteins) or a signal transduction ('non-classic'), proteins generally function involving stable or quasi-stable multi-protein associations. The physical basis for such associations is inherent in the combined effect of shape and electrostatic complementarities (Sc, EC) of the interacting protein partners at their interface, which provides indirect probabilistic estimates of the stability and affinity of the interaction. While Sc is a necessary criterion for inter-protein associations, EC can be favorable as well as disfavored (e.g., in transient interactions). Estimating equilibrium thermodynamic parameters (∆G_binding, K_d) by experimental means is costly and time consuming, thereby opening windows for computational structural interventions. Attempts to empirically probe ∆G_binding from coarse-grain structural descriptors (primarily, surface area based terms) have lately been overtaken by physics-based, knowledge-based and their hybrid approaches (MM/PBSA, FoldX, etc.) that directly compute ∆G_binding without involving intermediate structural descriptors.

METHODS

Here, we present EnCPdock ( https://www.scinetmol.in/EnCPdock/ ), a user-friendly web-interface for the direct conjoint comparative analyses of complementarity and binding energetics in proteins. EnCPdock returns an AI-predicted ∆G_binding computed by combining complementarity (Sc, EC) and other high-level structural descriptors (input feature vectors), and renders a prediction accuracy comparable to the state-of-the-art. EnCPdock further locates a PPI complex in terms of its {Sc, EC} values (taken as an ordered pair) in the two-dimensional complementarity plot (CP). In addition, it also generates mobile molecular graphics of the interfacial atomic contact network for further analyses. EnCPdock also furnishes individual feature trends along with the relative probability estimates (Pr_fmax) of the obtained feature-scores with respect to the events of their highest observed frequencies. Together, these functionalities are of real practical use for structural tinkering and intervention as might be relevant in the design of targeted protein-interfaces. Combining all its features and applications, EnCPdock presents a unique online tool that should be beneficial to structural biologists and researchers across related fraternities.

Additive energetic contributions of multiple peptide positions determine the relative promiscuity of viral and human sequences for PDZ domain targets

Protein-protein interactions that involve recognition of short peptides are critical in cellular processes. Protein-peptide interaction surface areas are relatively small and shallow, and there are often overlapping specificities in families of peptide-binding domains. Therefore, dissecting selectivity determinants can be challenging. PDZ domains are a family of peptide-binding domains located in several intracellular signaling and trafficking pathways. These domains are also directly targeted by pathogens, and a hallmark of many oncogenic viral proteins is a PDZ-binding motif. However, amidst sequences that target PDZ domains, there is a wide spectrum in relative promiscuity. For example, the viral HPV16 E6 oncoprotein recognizes over double the number of PDZ domain-containing proteins as the cystic fibrosis transmembrane conductance regulator (CFTR) in the cell, despite similar PDZ targeting-sequences and identical motif residues. Here, we determine binding affinities for PDZ domains known to bind either HPV16 E6 alone or both CFTR and HPV16 E6, using peptides matching WT and hybrid sequences. We also use energy minimization to model PDZ-peptide complexes and use sequence analyses to investigate this difference. We find that while the majority of single mutations had marginal effects on overall affinity, the additive effect on the free energy of binding accurately describes the selectivity observed. Taken together, our results describe how complex and differing PDZ interactomes can be programmed in the cell.

The Drug-Induced Interface That Drives HIV-1 Integrase Hypermultimerization and Loss of Function

Allosteric HIV-1 integrase (IN) inhibitors (ALLINIs) are an emerging class of small molecules that disrupt viral maturation by inducing the aberrant multimerization of IN. Here, we present cocrystal structures of HIV-1 IN with two potent ALLINIs, namely, BI-D and the drug candidate Pirmitegravir. The structures reveal atomistic details of the ALLINI-induced interface between the HIV-1 IN catalytic core and carboxyl-terminal domains (CCD and CTD). Projecting from their principal binding pocket on the IN CCD dimer, the compounds act as molecular glue by engaging a triad of invariant HIV-1 IN CTD residues, namely, Tyr226, Trp235, and Lys266, to nucleate the CTD-CCD interaction. The drug-induced interface involves the CTD SH3-like fold and extends to the beginning of the IN carboxyl-terminal tail region. We show that mutations of HIV-1 IN CTD residues that participate in the interface with the CCD greatly reduce the IN-aggregation properties of Pirmitegravir. Our results explain the mechanism of the ALLINI-induced condensation of HIV-1 IN and provide a reliable template for the rational development of this series of antiretrovirals through the optimization of their key contacts with the viral target. IMPORTANCE Despite the remarkable success of combination antiretroviral therapy, HIV-1 remains among the major causes of human suffering and loss of life in poor and developing nations. To prevail in this drawn-out battle with the pandemic, it is essential to continue developing advanced antiviral agents to fight drug resistant HIV-1 variants. Allosteric integrase inhibitors (ALLINIs) are an emerging class of HIV-1 antagonists that are orthogonal to the current antiretroviral drugs. These small molecules act as highly specific molecular glue, which triggers the aggregation of HIV-1 integrase. In this work, we present high-resolution crystal structures that reveal the crucial interactions made by two potent ALLINIs, namely, BI-D and Pirmitegravir, with HIV-1 integrase. Our results explain the mechanism of drug action and will inform the development of this promising class of small molecules for future use in antiretroviral regimens.

Alchemical Free Energy Calculations to Investigate Protein-Protein Interactions: the Case of the CDC42/PAK1 Complex

Here, we show that alchemical free energy calculations can quantitatively compute the effect of mutations at the protein–protein interface. As a test case, we have used the protein complex formed by the small Rho-GTPase CDC42 and its downstream effector PAK1, a serine/threonine kinase. Notably, the CDC42/PAK1 complex offers a wealth of structural, mutagenesis, and binding affinity data because of its central role in cellular signaling and cancer progression. In this context, we have considered 16 mutations in the CDC42/PAK1 complex and obtained excellent agreement between computed and experimental data on binding affinity. Importantly, we also show that a careful analysis of the side-chain conformations in the mutated amino acids can considerably improve the computed estimates, solving issues related to sampling limitations. Overall, this study demonstrates that alchemical free energy calculations can conveniently be integrated into the design of experimental mutagenesis studies.

How Peptides Bind to PSD-95/Discs-Large/ZO-1 Domains

PSD-95/discs-large/ZO-1 (PDZ) domains form a large family of adaptor proteins that bind to the C-terminal tails of their binding partner proteins. Via extensive molecular dynamics simulations and alchemical free energy calculations, we characterized the binding modi of phosphorylated and unphosphorylated EQVSAV peptides and of a EQVEAV phosphate mimic to the hPTP1E PDZ2 and MAGI1 PDZ1 domains. The simulations reproduced the well-known binding characteristics such as tight coordination of the peptidic carboxyl tail and pronounced hydrogen bonding between the peptide backbone and the backbone atoms of a β-sheet in PDZ. Overall, coordination by hPTP1E PDZ2 appeared tighter than by MAGI1 PDZ1. Simulations of wild-type PDZ and arginine mutants suggest that contacts with Arg79/85 in hPTP1E/MAGI1 are more important for the EQVEAV peptide than for EQVSAV. Alchemical free energy calculations and PaCS-MD simulations could well reproduce the difference in binding free energy between unphosphorylated EQVSAV and EQVEAV peptides and the absolute binding free energy of EQVSAV. However, likely due to small force field inaccuracies, the simulations erroneously favored binding of the phosphorylated peptide instead of its unphosphorylated counterpart, which is in contrast to the experiment.

Dawn of a New Era for Membrane Protein Design

A major advancement has recently occurred in the ability to predict protein secondary structure from sequence using artificial neural networks. This new accessibility to high-quality predicted structures provides a big opportunity for the protein design community. It is particularly welcome for membrane protein design, where the scarcity of solved structures has been a major limitation of the field for decades. Here, we review the work done to date on the membrane protein design and set out established and emerging tools that can be used to most effectively exploit this new access to structures.

Integrated unsupervised-supervised modeling and prediction of protein-peptide affinities at structural level

Cell signal networks are orchestrated directly or indirectly by various peptide-mediated protein-protein interactions, which are normally weak and transient and thus ideal for biological regulation and medicinal intervention. Here, we develop a general-purpose method for modeling and predicting the binding affinities of protein-peptide interactions (PpIs) at the structural level. The method is a hybrid strategy that employs an unsupervised approach to derive a layered PpI atom-residue interaction (ulPpI[a-r]) potential between different protein atom types and peptide residue types from thousands of solved PpI complex structures and then statistically correlates the potential descriptors with experimental affinities (KD values) over hundreds of known PpI samples in a supervised manner to create an integrated unsupervised-supervised PpI affinity (usPpIA) predictor. Although both the ulPpI[a-r] potential and usPpIA predictor can be used to calculate PpI affinities from their complex structures, the latter seems to perform much better than the former, suggesting that the unsupervised potential can be improved substantially with a further correction by supervised statistical learning. We examine the robustness and fault-tolerance of usPpIA predictor when applied to treat the coarse-grained PpI complex structures modeled computationally by sophisticated peptide docking and dynamics simulation. It is revealed that, despite developed solely based on solved structures, the integrated unsupervised-supervised method is also applicable for locally docked structures to reach a quantitative prediction but can only give a qualitative prediction on globally docked structures. The dynamics refinement seems not to change (or improve) the predictive results essentially, although it is computationally expensive and time-consuming relative to peptide docking. We also perform extrapolation of usPpIA predictor to the indirect affinity quantities of HLA-A*0201 binding epitope peptides and NHERF PDZ binding scaffold peptides, consequently resulting in a good and moderate correlation of the predicted KD with experimental IC50 and BLU on the two peptide sets, with Pearson's correlation coefficients Rp = 0.635 and 0.406, respectively.

Relative Binding Free Energy Calculations for Ligands with Diverse Scaffolds with the Alchemical Transfer Method

We present an extension of the alchemical transfer method (ATM) for the estimation of relative binding free energies of molecular complexes applicable to conventional, as well as scaffold-hopping, alchemical transformations. Named ATM-RBFE, the method is implemented in the free and open-source OpenMM molecular simulation package and aims to provide a simpler and more generally applicable route to the calculation of relative binding free energies than what is currently available. ATM-RBFE is based on sound statistical mechanics theory and a novel coordinate perturbation scheme designed to swap the positions of a pair of ligands such that one is transferred from the bulk solvent to the receptor binding site while the other moves simultaneously in the opposite direction. The calculation is conducted directly in a single solvent box with a system prepared with conventional setup tools, without splitting of electrostatic and nonelectrostatic transformations, and without pairwise soft-core potentials. ATM-RBFE is validated here against the absolute binding free energies of the SAMPL8 GDCC host-guest benchmark set and against protein-ligand benchmark sets that include complexes of the estrogen receptor ERα and those of the methyltransferase EZH2. In each case the method yields self-consistent and converged relative binding free energy estimates in agreement with absolute binding free energies and reference literature values, as well as experimental measurements.

Free Energy-Based Computational Methods for the Study of Protein-Peptide Binding Equilibria

This chapter discusses the theory and application of physics-based free energy methods to estimate protein-peptide binding free energies. It presents a statistical mechanics formulation of molecular binding, which is then specialized in three methodologies: (1) alchemical absolute binding free energy estimation with implicit solvation, (2) alchemical relative binding free energy estimation with explicit solvation, and (3) potential of mean force binding free energy estimation. Case studies of protein-peptide binding application taken from the recent literature are discussed for each method.

Fast Calculation of Electrostatic Solvation Free Energy in Simple Ionic Fluids Using an Energy-Scaled Debye-Hückel Theory

Continuum theories are useful to compute the solvation free energy in ionic fluids. Herein, the electrostatic solvation free energy (ESFE) in simple ionic fluids is studied with an energy-scaled Debye-Hückel (ESDH) theory. Given the ESFEs of simple spherical ions as input, the ESDH theory is applicable to molecules with various complex geometries and charge distributions. Specifically, the ESDH theory is applied to molecules in a molten salt system, where the predicted ESFEs are in good agreement with molecular dynamics simulation results. Our study sheds light on accurately predicting the ESFE in ionic fluids with phenomenological continuum theories.

Alchemical Transfer Approach to Absolute Binding Free Energy Estimation

The alchemical transfer method (ATM) for the calculation of standard binding free energies of noncovalent molecular complexes is presented. The method is based on a coordinate displacement perturbation of the ligand between the receptor binding site and the explicit solvent bulk and a thermodynamic cycle connected by a symmetric intermediate in which the ligand interacts with the receptor and solvent environments with equal strength. While the approach is alchemical, the implementation of the ATM is as straightforward as that for physical pathway methods of binding. The method is applicable, in principle, with any force field, as it does not require splitting the alchemical transformations into electrostatic and nonelectrostatic steps, and it does not require soft-core pair potentials. We have implemented the ATM as a freely available and open-source plugin of the OpenMM molecular dynamics library. The method and its implementation are validated on the SAMPL6 SAMPLing host-guest benchmark set. The work paves the way to streamlined alchemical relative and absolute binding free energy implementations on many molecular simulation packages and with arbitrary energy functions including polarizable, quantum-mechanical, and artificial neural network potentials.

Computational Design of PDZ-Peptide Binding

This chapter describes two computational methods for PDZ-peptide binding: high-throughput computational protein design (CPD) and a medium-throughput approach combining molecular dynamics for conformational sampling with a Poisson-Boltzmann (PB) Linear Interaction Energy for scoring. A new CPD method is outlined, which uses adaptive Monte Carlo simulations to efficiently sample peptide variants that tightly bind a PDZ domain, and provides at the same time precise estimates of their relative binding free energies. A detailed protocol is described based on the Proteus CPD software. The medium-throughput approach can be performed with standard MD and PB software, such as NAMD and Charmm. For 40 complexes between Tiam1 and peptide ligands, it gave high a2ccuracy, with mean errors of around 0.5 kcal/mol for relative binding free energies and no large errors. It requires a moderate amount of parameter fitting before it can be applied, and its transferability to other protein families is still untested.

Statistical mechanics of polarizable force fields based on classical Drude oscillators with dynamical propagation by the dual-thermostat extended Lagrangian

Polarizable force fields based on classical Drude oscillators offer a practical and computationally efficient avenue to carry out molecular dynamics (MD) simulations of large biomolecular systems. To treat the polarizable electronic degrees of freedom, the Drude model introduces a virtual charged particle that is attached to its parent nucleus via a harmonic spring. Traditionally, the need to relax the electronic degrees of freedom for each fixed set of nuclear coordinates is achieved by performing an iterative self-consistent field (SCF) calculation to satisfy a selected tolerance. This is a computationally demanding procedure that can increase the computational cost of MD simulations by nearly one order of magnitude. To avoid the costly SCF procedure, a small mass is assigned to the Drude particles, which are then propagated as dynamic variables during the simulations via a dual-thermostat extended Lagrangian algorithm. To help clarify the significance of the dual-thermostat extended Lagrangian propagation in the context of the polarizable force field based on classical Drude oscillators, the statistical mechanics of a dual-temperature canonical ensemble is formulated. The conditions for dynamically maintaining the dual-temperature properties in the case of the classical Drude oscillator are analyzed using the generalized Langevin equation.

Improved Modeling of Cation-π and Anion-Ring Interactions Using the Drude Polarizable Empirical Force Field for Proteins

Cation-π interactions are noncovalent interactions between a π-electron system and a positively charged ion that are regarded as a strong noncovalent interaction and are ubiquitous in biological systems. Similarly, though less studied, anion-ring interactions are present in proteins along with in-plane interactions of anions with aromatic rings. As these interactions are between a polarizing ion and a polarizable π system, the accuracy of the treatment of these interactions in molecular dynamics (MD) simulations using additive force fields (FFs) may be limited. In the present work, to allow for a better description of ion-π interactions in proteins in the Drude-2013 protein polarizable FF, we systematically optimized the parameters for these interactions targeting model compound quantum mechanical (QM) interaction energies with atom pair-specific Lennard-Jones parameters along with virtual particles as selected ring centroids introduced to target the QM interaction energies and geometries. Subsequently, MD simulations were performed on a series of protein structures where ion-π pairs occur to evaluate the optimized parameters in the context of the Drude-2013 FF. The resulting FF leads to a significant improvement in reproducing the ion-π pair distances observed in experimental protein structures, as well as a smaller root-mean-square differences and fluctuations of the overall protein structures from experimental structures. Accordingly, the optimized Drude-2013 protein polarizable FF is suggested for use in MD simulations of proteins where cation-π and anion-ring interactions are critical. © 2019 Wiley Periodicals, Inc.

Polarisable force fields: what do they add in biomolecular simulations?

The quality of biomolecular simulations critically depends on the accuracy of the force field used to calculate the potential energy of the molecular configurations. Currently, most simulations employ non-polarisable force fields, which describe electrostatic interactions as the sum of Coulombic interactions between fixed atomic charges. Polarisation of these charge distributions is incorporated only in a mean-field manner. In the past decade, extensive efforts have been devoted to developing simple, efficient, and yet generally applicable polarisable force fields for biomolecular simulations. In this review, we summarise the latest developments in accounting for key biomolecular interactions with polarisable force fields and applications to address challenging biological questions. In the end, we provide an outlook for future development in polarisable force fields.

Rapid and accurate structure-based therapeutic peptide design using GPU accelerated thermodynamic integration

Peptide-based therapeutics are an alternative to small molecule drugs as they offer superior specificity, lower toxicity, and easy synthesis. Here we present an approach that leverages the dramatic performance increase afforded by the recent arrival of GPU accelerated thermodynamic integration (TI). GPU TI facilitates very fast, highly accurate binding affinity optimization of peptides against therapeutic targets. We benchmarked TI predictions using published peptide binding optimization studies. Prediction of mutations involving charged side-chains was found to be less accurate than for non-charged, and use of a more complex 3-step TI protocol was found to boost accuracy in these cases. Using the 3-step protocol for non-charged side-chains either had no effect or was detrimental. We use the benchmarked pipeline to optimize a peptide binding to our recently discovered cancer target: EME1. TI calculations predict beneficial mutations using both canonical and non-canonical amino acids. We validate these predictions using fluorescence polarization and confirm that binding affinity is increased. We further demonstrate that this increase translates to a significant reduction in pancreatic cancer cell viability.

Adaptive landscape flattening in amino acid sequence space for the computational design of protein:peptide binding

For the high throughput design of protein:peptide binding, one must explore a vast space of amino acid sequences in search of low binding free energies. This complex problem is usually addressed with either simple heuristic scoring or expensive sequence enumeration schemes. Far more efficient than enumeration is a recent Monte Carlo approach that adaptively flattens the energy landscape in sequence space of the unbound peptide and provides formally exact binding free energy differences. The method allows the binding free energy to be used directly as the design criterion. We propose several improvements that allow still more efficient sampling and can address larger design problems. They include the use of Replica Exchange Monte Carlo and landscape flattening for both the unbound and bound peptides. We used the method to design peptides that bind to the PDZ domain of the Tiam1 signaling protein and could serve as inhibitors of its activity. Four peptide positions were allowed to mutate freely. Almost 75 000 peptide variants were processed in two simulations of 10⁹ steps each that used 1 CPU hour on a desktop machine. 96% of the theoretical sequence space was sampled. The relative binding free energies agreed qualitatively with values from experiment. The sampled sequences agreed qualitatively with an experimental library of Tiam1-binding peptides. The main assumption limiting accuracy is the fixed backbone approximation, which could be alleviated in future work by using increased computational resources and multi-backbone designs.

Classical Drude Polarizable Force Field Model for Methyl Phosphate and Its Interactions with Mg2

Phosphate groups are essential components of nucleic acids and proteins, whose interactions with solvent, metal ions, and ionic side chains help control folding and binding. Methyl phosphate (MP) represents a simple analog of phosphate moieties that are post-translation modifications in proteins and present at the termini of nucleic acids, among other environments. In the present study, we optimized parameters for use in polarizable molecular dynamics simulations of MP in its mono- and dianionic forms, MP- ≡ CH3HPO4- and MP2- ≡ CH3PO42-, along with P i2- ≡ HPO42-, in the context of the classical Drude oscillator model. Parameter optimization was done in a manner consistent with the remainder of the Drude molecular mechanics force field, choosing atomic charges and polarizabilities to reproduce molecular properties from quantum mechanics as well as experimental hydration free energies. Optimized parameters were similar to existing dimethyl phosphate parameters, with a few significant differences. The developed parameters were then used to compute magnesium binding affinities in aqueous solution, using alchemical molecular dynamics free energy simulations. Good agreement with experiment was obtained, and outer sphere binding was shown to be predominant for MP- and MP2-.