Relative contributions of desolvation, inter- and intramolecular interactions to binding affinity in protein kinase systems

In several previous studies, we performed sensitivity analysis to gauge the relative importance of different atomic partial charges in determining protein-ligand binding. In this work, we gain further insights by decomposing these results into three contributions: desolvation, intramolecular interactions, and intermolecular interactions, again based on a Poisson continuum electrostatics model. Three protein kinase-inhibitor systems have been analyzed: CDK2-deschloroflavopiridol, PKA-PKI, and LCK-PP2. Although our results point out the importance of specific intermolecular interactions to the binding affinity, they also reveal the remarkable contributions from the solvent-mediated intramolecular interactions in some cases. Thus, it is necessary to look beyond analyzing protein-ligand interactions to understand protein-ligand recognition or to gain insights into designing ligands and proteins. In analyzing the contributions of the three components to the overall binding free energy, the PKA-PKI system with a much larger ligand was found to behave differently from the other two systems with smaller ligands. In the former case, the intermolecular interactions are very favorable, and together with the favorable solvent-mediated intramolecular interactions, they overcome the large desolvation penalties to give a favorable electrostatics contribution to the overall binding affinity. On the other hand, the other two systems with smaller ligands only present modest intermolecular interactions and they are not or are only barely sufficient to overcome the desolvation penalty even with the aid of the favorable intramolecular contributions. As a result, the binding affinity of these two systems do not or only barely benefit from electrostatics contributions.     

Calculation of ligand binding free energies from molecular dynamics simulations

A recently developed method for predicting binding affinities in ligand–receptor complexes, based on interaction energy averaging and conformational sampling by molecular dynamics simulation, is presented. Polar and nonpolar contributions to the binding free energy are approximated by a linear scaling of the corresponding terms in the average intermolecular interaction energy for the bound and free states of the ligand. While the method originally assumed the validity of electrostatic linear response, we show that incorporation of systematic deviations from linear response derived from free energy perturbation calculations enhances the accuracy of the approach. The method is applied to complexes of wild-type and mutant human dihydrofolate reductases with 2,4-diaminopteridine and 2,4-diaminoquinazoline inhibitors. It is shown that a binding energy accuracy of about 1 kcal/mol is attainable even for multiply ionized compounds, such as methotrexate, for which electrostatic interactions energies are very large. © 1998 John Wiley & Sons, Inc. Int J Quant Chem 69: 77–88, 1998     

Influencing receptor-ligand binding mechanisms with multivalent ligand architecture

Multivalent ligands can function as inhibitors or effectors of biological processes. Potent inhibitory activity can arise from the high functional affinities of multivalent ligand-receptor interactions. Effector functions, however, are influenced not only by apparent affinities but also by alternate factors, including the ability of a ligand to cluster receptors. Little is known about the molecular features of a multivalent ligand that determine whether it will function as an inhibitor or effector. We envisioned that, by altering multivalent ligand architecture, ligands with preferences for different binding mechanisms would be generated. To this end, a series of 28 ligands possessing structural diversity was synthesized. This series provides the means to explore the effects of ligand architecture on the inhibition and clustering of a model protein, the lectin concanavalin A (Con A). The structural parameters that were varied include scaffold shape, size, valency, and density of binding elements. We found that ligands with certain architectures are effective inhibitors, but others mediate receptor clustering. Specifically, high molecular weight, polydisperse polyvalent ligands are effective inhibitors of Con A binding, whereas linear oligomeric ligands generated by the ring-opening metathesis polymerization have structural properties that favor clustering. The shape of a multivalent ligand also influences specific aspects of receptor clustering. These include the rate at which the receptor is clustered, the number of receptors in the clusters, and the average interreceptor distance. Our results indicate that the architecture of a multivalent ligand is a key parameter in determining its activity as an inhibitor or effector. Diversity-oriented syntheses of multivalent ligands coupled with effective assays that can be used to compare the contributions of different binding parameters may afford ligands that function by specific mechanisms.     

Modeling multiple species of nicotine and deschloroepibatidine interacting with alpha4beta2 nicotinic acetylcholine receptor: from microscopic binding to phenomenological binding affinity

A variety of molecular modeling, molecular docking, and first-principles electronic structure calculations were performed to study how the alpha4beta2 nicotinic acetylcholine receptor (nAChR) binds with different species of two typical agonists, (S)-(-)-nicotine and (R)-(-)-deschloroepibatidine, each of which is distinguished by different free bases and protonation states. On the basis of these results, predictions were made regarding the corresponding microscopic binding free energies. Hydrogen-bonding and cation-pi interactions between the receptor and the respective ligands were found to be the dominant factors differentiating the binding strengths of different microscopic binding species. The calculated results and analyses demonstrate that, for each agonist, all the species are interchangeable and can quickly achieve a thermodynamic equilibrium in solution and at the nAChR binding site. This allows quantitation of the equilibrium concentration distributions of the free ligand species and the corresponding microscopic ligand-receptor binding species, their pH dependence, and their contributions to the phenomenological binding affinity. The predicted equilibrium concentration distributions, pK(a) values, absolute phenomenological binding affinities, and their pH dependence are all in good agreement with available experimental data, suggesting that the computational strategy from the microscopic binding species and affinities to the phenomenological binding affinity is reliable for studying alpha4beta2 nAChR-ligand binding. This should provide valuable information for future rational design of drugs targeting nAChRs. The general strategy of the "from-microscopic-to-phenomenological" approach for studying interactions of alpha4beta2 nAChRs with (S)-(-)-nicotine and (R)-(-)-deschloroepibatidine may also be useful in studying other types of ligand-protein interactions involving multiple molecular species of a ligand and in associated rational drug design.     

Thermodynamics of buried water clusters at a protein-ligand binding interface

The structure of the complex of cyclophilin A (CypA) with cyclosporin A (CsA, 1) shows a cluster of four water molecules buried at the binding interface, which is rearranged when CsA is replaced by (5-hydroxynorvaline)-2-cyclosporin (2). The thermodynamic contributions of each bound water molecule in the two complexes are explored with the inhomogeneous fluid solvation theory and molecular dynamics simulations. Water (WTR) 133 in complex 1 contributes little to the binding affinity, while WTR6 and 7 in complex 2 play an essential role in mediating protein-ligand binding with a hydrogen bond network. The calculations reveal that the rearrangement of the water molecules contributes favorably to the binding affinity, even though one of them is displaced going from ligand 1 to 2. Another favorable contribution comes from the larger protein-ligand interactions of ligand 2. However, these favorable contributions are not sufficient to overcome the unfavorable desolvation free energy change and the conformational entropy of the hydroxylpropyl group of ligand 2 in the complex, leading to a lower binding affinity of ligand 2. These physical insights may be useful in the development of improved scoring functions for binding affinity prediction.     

Prediction of trypsin/molecular fragment binding affinities by free energy decomposition and empirical scores

Two families of binding affinity estimation methodologies are described which were utilized in the SAMPL3 trypsin/fragment binding affinity challenge. The first is a free energy decomposition scheme based on a thermodynamic cycle, which included separate contributions from enthalpy and entropy of binding as well as a solvent contribution. Enthalpic contributions were estimated with PM6-DH2 semiempirical quantum mechanical interaction energies, which were modified with a statistical error correction procedure. Entropic contributions were estimated with the rigid-rotor harmonic approximation, and solvent contributions to the free energy were estimated with several different methods. The second general methodology is the empirical score LISA, which contains several physics-based terms trained with the large PDBBind database of protein/ligand complexes. Here we also introduce LISA+, an updated version of LISA which, prior to scoring, classifies systems into one of four classes based on a ligand's hydrophobicity and molecular weight. Each version of the two methodologies (a total of 11 methods) was trained against a compiled set of known trypsin binders available in the Protein Data Bank to yield scaling parameters for linear regression models. Both raw and scaled scores were submitted to SAMPL3. Variants of LISA showed relatively low absolute errors but also low correlation with experiment, while the free energy decomposition methods had modest success when scaling factors were included. Nonetheless, re-scaled LISA yielded the best predictions in the challenge in terms of RMS error, and six of these models placed in the top ten best predictions by RMS error. This work highlights some of the difficulties of predicting binding affinities of small molecular fragments to protein receptors as well as the benefit of using training data.     

Sensitivity to antitumor drugs and vinblastine binding to membrane in rat ascites hepatoma AH66 cells.

Rat ascites hepatoma AH66 cells have lower sensitivity to Vinca alkaloids and anthracycline antibiotics than AH66F cells, a subline of AH66 cells. AH66 cells expressed P-glycoprotein, while the protein was not detectable in AH66F cells. There are two affinity sites for [3H]vinblastine binding in the AH66 cell membrane, while AH66F cells have only one affinity site. The high affinity [3H]vinblastine binding in AH66 cells was inhibited by Adriamycin, verapamil, nicardipine, and reserpine. The high affinity site of the binding may be the multidrug transporter, P-glycoprotein. [3H]Vinblastine binding was not influenced by adenosine 3'-5'-monophosphate (AMP), adenosine triphosphate (ATP), or guanosine triphosphate (GTP). The multidrug resistance in AH66 cells may depend on P-glycoprotein which is not modulated by nucleotide.     

Rapid and accurate prediction of binding free energies for saquinavir-bound HIV-1 proteases

To explain drug resistance by computer simulations at the molecular level, we first have to assess the accuracy of theoretical predictions. Herein we report an application of the molecular mechanics Poisson-Boltzmann surface area (MM/PBSA) technique to the ranking of binding affinities of the inhibitor saquinavir with the wild type (WT) and three resistant mutants of HIV-1 protease: L90M, G48V, and G48V/L90M. For each ligand-protein complex we report 10 ns of fully unrestrained molecular dynamics (MD) simulations with explicit solvent. We investigate convergence, internal consistency, and model dependency of MM/PBSA ligand binding energies. Converged enthalpy and entropy estimates produce ligand binding affinities within 1.5 kcal/mol of experimental values, with a remarkable level of correlation to the experimentally observed ranking of resistance levels. A detailed analysis of the enthalpic/entropic balance of drug-protease interactions explains resistance in L90M in terms of a higher vibrational entropy than in the WT complex, while G48V disrupts critical hydrogen bonds at the inhibitor's binding site and produces an altered, more unfavorable balance of Coulomb and polar desolvation energies.     

Processing multimode binding situations in simulation-based prediction of ligand-macromolecule affinities

The linear response (LR) approximation and similar approaches belong to practical methods for estimation of ligand-receptor binding affinities. The approaches correlate experimental binding affinities with the changes upon binding of the ligand electrostatic and van der Waals energies and of solvation characteristics. These attributes are expressed as ensemble averages that are obtained by conformational sampling of the protein-ligand complex and of the free ligand by molecular dynamics or Monte Carlo simulations. We observed that outliers in the LR correlations occasionally exhibit major conformational changes of the complex during sampling. We treated the situation as a multimode binding case, for which the observed association constant is the sum of the partial association constants of individual states/modes. The resulting nonlinear expression for the binding affinities contains all the LR variables for individual modes that are scaled by the same two to four adjustable parameters as in the one-mode LR equation. The multimode method was applied to inhibitors of a matrix metalloproteinase, where this treatment improved the explained variance in experimental activity from 75% for the unimode case to about 85%. The predictive ability scaled accordingly, as verified by extensive cross-validations.     

Dual Stimuli-Responsive Controlled Release Nanocarrier for Multidrug Resistance Cancer Therapy

Multidrug resistance of cancer cells is a major obstacle for cancer chemotherapy. Herein, we present a nanocarrier that can release chemotherapeutic agents to induce tumor cell death and generate NO under NIR to overcome multidrug resistance in cancer chemotherapy. Owing to the unique structure of the water channel in this controlled release system for chemotherapeutic agents, the nanocarrier surface is equipped with more active sites to graft NO donor molecules. The released NO performs very well in reversing multidrug resistance by inhibiting P-gp expression. Our findings provide new insight into multidrug resistance cancer therapy and controlled release nanocarriers for multiple drugs.     

Rapid and accurate estimation of protein–ligand relative binding affinities using site-identification by ligand competitive saturation

Predicting relative protein–ligand binding affinities is a central pillar of lead optimization efforts in structure-based drug design. The site identification by ligand competitive saturation (SILCS) methodology is based on functional group affinity patterns in the form of free energy maps that may be used to compute protein–ligand binding poses and affinities. Presented are results obtained from the SILCS methodology for a set of eight target proteins as reported originally in Wang et al. (J. Am. Chem. Soc., 2015, 137, 2695–2703) using free energy perturbation (FEP) methods in conjunction with enhanced sampling and cycle closure corrections. These eight targets have been subsequently studied by many other authors to compare the efficacy of their method while comparing with the outcomes of Wang et al. In this work, we present results for a total of 407 ligands on the eight targets and include specific analysis on the subset of 199 ligands considered previously. Using the SILCS methodology we can achieve an average accuracy of up to 77% and 74% when considering the eight targets with their 199 and 407 ligands, respectively, for rank-ordering ligand affinities as calculated by the percent correct metric. This accuracy increases to 82% and 80%, respectively, when the SILCS atomic free energy contributions are optimized using a Bayesian Markov-chain Monte Carlo approach. We also report other metrics including Pearson''s correlation coefficient, Pearlman''s predictive index, mean unsigned error, and root mean square error for both sets of ligands. The results obtained for the 199 ligands are compared with the outcomes of Wang et al. and other published works. Overall, the SILCS methodology yields similar or better-quality predictions without a priori need for known ligand orientations in terms of the different metrics when compared to current FEP approaches with significant computational savings while additionally offering quantitative estimates of individual atomic contributions to binding free energies. These results further validate the SILCS methodology as an accurate, computationally efficient tool to support lead optimization and drug discovery.

Predicting relative protein–ligand binding affinities is a central pillar of lead optimization efforts in structure-based drug design.     

Extended-ensemble docking to probe dynamic variation of ligand binding sites during large-scale structural changes of proteins

Proteins can sample a broad landscape as they undergo conformational transition between different functional states. At the same time, as key players in almost all cellular processes, proteins are important drug targets. Considering the different conformational states of a protein is therefore central for a successful drug-design strategy. Here we introduce a novel docking protocol, termed extended-ensemble docking, pertaining to proteins that undergo large-scale (global) conformational changes during their function. In its application to multidrug ABC-transporter P-glycoprotein (Pgp), extensive non-equilibrium molecular dynamics simulations employing system-specific collective variables are first used to describe the transition cycle of the transporter. An extended set of conformations (extended ensemble) representing the full transition cycle between the inward- and the outward-facing states is then used to seed high-throughput docking calculations of known substrates, non-substrates, and modulators of the transporter. Large differences are predicted in the binding affinities to different conformations, with compounds showing stronger binding affinities to intermediate conformations compared to the starting crystal structure. Hierarchical clustering of the binding modes shows all ligands preferably bind to the large central cavity of the protein, formed at the apex of the transmembrane domain (TMD), whereas only small binding populations are observed in the previously described R and H sites present within the individual TMD leaflets. Based on the results, the central cavity is further divided into two major subsites, first preferably binding smaller substrates and high-affinity inhibitors, whereas the second one shows preference for larger substrates and low-affinity modulators. These central subsites along with the low-affinity interaction sites present within the individual TMD leaflets may respectively correspond to the proposed high- and low-affinity binding sites in Pgp. We propose further an optimization strategy for developing more potent inhibitors of Pgp, based on increasing its specificity to the extended ensemble of the protein, instead of using a single protein structure, as well as its selectivity for the high-affinity binding site. In contrast to earlier in silico studies using single static structures of Pgp, our results show better agreement with experimental studies, pointing to the importance of incorporating the global conformational flexibility of proteins in future drug-discovery endeavors.

Functional states of P-glycoprotein formed during its full transition cycle (red to blue), captured by molecular dynamics simulations, form a structural framework for extended-ensemble docking of small-molecule ligands of diverse activities.     

The Structural Basis of Mycobacterium tuberculosis RpoB Drug-Resistant Clinical Mutations on Rifampicin Drug Binding

Tuberculosis (TB), caused by the Mycobacterium tuberculosis infection, continues to be a leading cause of morbidity and mortality in developing countries. Resistance to the first-line anti-TB drugs, isoniazid (INH) and rifampicin (RIF), is a major drawback to effective TB treatment. Genetic mutations in the β-subunit of the DNA-directed RNA polymerase (rpoB) are reported to be a major reason of RIF resistance. However, the structural basis and mechanisms of these resistant mutations are insufficiently understood. In the present study, thirty drug-resistant mutants of rpoB were initially modeled and screened against RIF via a comparative molecular docking analysis with the wild-type (WT) model. These analyses prioritized six mutants (Asp441Val, Ser456Trp, Ser456Gln, Arg454Gln, His451Gly, and His451Pro) that showed adverse binding affinities, molecular interactions, and RIF binding hinderance properties, with respect to the WT. These mutant models were subsequently analyzed by molecular dynamics (MD) simulations. One-hundred nanosecond all-atom MD simulations, binding free energy calculations, and a dynamic residue network analysis (DRN) were employed to exhaustively assess the impact of mutations on RIF binding dynamics. Considering the global structural motions and protein–ligand binding affinities, the Asp441Val, Ser456Gln, and His454Pro mutations generally yielded detrimental effects on RIF binding. Locally, we found that the electrostatic contributions to binding, particularly by Arg454 and Glu487, might be adjusted to counteract resistance. The DRN analysis revealed that all mutations mostly distorted the communication values of the critical hubs and may, therefore, confer conformational changes in rpoB to perturb RIF binding. In principle, the approach combined fundamental molecular modeling tools for robust “global” and “local” level analyses of structural dynamics, making it well suited for investigating other similar drug resistance cases.     

Computational study of fluoroquinolone binding to and its applicability to future drug design

Fluoroquinolones are an important therapeutic class in the targeting of new and resistant bacterial infections. Fluoroquinolones bind to bacterial type II topoisomerase via a water‐Mg²⁺ bridge. However, binding to magnesium‐containing molecules outside of the target cells increases the minimum inhibitory concentration (MIC) and promotes drug resistance. As a result, fluoroquinolones are counter‐indicated with magnesium and multivalent metal cation containing drugs, such as antacids. The antibiotic efficacy of fluoroquinolones has also been shown to be pH dependent, as we show the effect of protonation state on magnesium binding. This work presents a systematic computational study of fluoroquinolones' magnesium‐binding properties. We use B3LYP density functional theory and triple‐zeta basis sets, to evaluate binding affinities. Complexation is predicted to be thermodynamically favorable at neutral and basic compared to acidic pH. The calculated complexation energies broadly capture experimental binding affinities, suggesting this is a valid approach for designing new fluoroquinolones with a target magnesium binding affinity. We also investigate the effect of chemical substitution at the carboxylic acid to help in the identification of potential new antibiotics based on the fluoroquinolone pharmacophore.     

Supramolecular binding thermodynamics by dispersion-corrected density functional theory

The equilibrium association free enthalpies ΔG(a) for typical supramolecular complexes in solution are calculated by ab initio quantum chemical methods. Ten neutral and three positively charged complexes with experimental ΔG(a) values in the range 0 to -21?kcal?mol(-1) (on average -6?kcal?mol(-1) ) are investigated. The theoretical approach employs a (nondynamic) single-structure model, but computes the various energy terms accurately without any special empirical adjustments. Dispersion corrected density functional theory (DFT-D3) with extended basis sets (triple-ζ and quadruple-ζ quality) is used to determine structures and gas-phase interaction energies (ΔE), the COSMO-RS continuum solvation model (based on DFT data) provides solvation free enthalpies and the remaining ro-vibrational enthalpic/entropic contributions are obtained from harmonic frequency calculations. Low-lying vibrational modes are treated by a free-rotor approximation. The accurate account of London dispersion interactions is mandatory with contributions in the range -5 to -60?kcal?mol(-1) (up to 200?% of ΔE). Inclusion of three-body dispersion effects improves the results considerably. A semilocal (TPSS) and a hybrid density functional (PW6B95) have been tested. Although the ΔG(a) values result as a sum of individually large terms with opposite sign (ΔE vs. solvation and entropy change), the approach provides unprecedented accuracy for ΔG(a) values with errors of only 2?kcal?mol(-1) on average. Relative affinities for different guests inside the same host are always obtained correctly. The procedure is suggested as a predictive tool in supramolecular chemistry and can be applied routinely to semirigid systems with 300-400 atoms. The various contributions to binding and enthalpy-entropy compensations are discussed.     

Amide- and urea-functionalized dithienylethene: synthesis, photochromism, and binding with halide anions

A versatile amide- and urea-functionalized dithienylethene has been successfully synthesized. Upon irradiation with UV or visible light, the compound showed excellent fatigue resistance. As a synthetic receptor, the dithienylethene displayed switchable affinities for Cl(-) and Br(-) anions when the UV/vis light was introduced. The switchable binding ability also had good reversibility.     

Water-soluble amino derivatives of free-base dppz--syntheses and DNA binding studies

The syntheses of two water soluble dipyrido-[3,2-a:2',3'-c]-phenazine analogues containing one or two appended amino/amide chains are reported. Steady state optical studies on the two new compounds reveal high-energy dppz-based luminescence in water and non-aqueous solvents. Optical titrations with duplex DNA show that the luminescence is quenched on the addition of DNA. Binding curves constructed from absorption and emission changes indicate that, while one of the compounds display negligible binding properties, the other binds DNA with relatively high affinity (>10(5) M(-1)). Isothermal calorimetry experiments, designed to investigate the higher binding compound in more detail, reveals that its interaction with CT-DNA is actually biphasic with one tight (>10(5) M(-1)) and one weaker binding site (~10(5) M(-1)). In both cases binding is entropically driven. Further calorimetry studies involving the interaction of the new compound with a variety of polynucleotides were carried out. To aid comparisons, similar experiments involving a previously reported bipyridyldiylium derivative of dppz were also carried out. These studies reveal that the bipyridyldiylium derivative binds all these sequences monophasically with relatively low affinities (~10(4) M(-1)). However, while the amino/amide chain appended derivative binds to Poly(dA).poly(dT) monophasically with relatively low affinities, it binds all the other polynucleotide studied biphasically, with affinities ranging from <10(6) M(-1) to >10(8) M(-1). The ITC data reveals that for both compounds thermodynamic signatures for binding are dependent on the sequence being bound. In both cases, the data for Poly(dA).poly(dT) is particularly anomalous. An analysis of the data shows that binding is selective, with affinities at flexible sequences being several orders of magnitude higher than those at more rigid sequences.