The free energy landscape of small molecule unbinding

The spontaneous dissociation of six small ligands from the active site of FKBP (the FK506 binding protein) is investigated by explicit water molecular dynamics simulations and network analysis. The ligands have between four (dimethylsulphoxide) and eleven (5-diethylamino-2-pentanone) non-hydrogen atoms, and an affinity for FKBP ranging from 20 to 0.2 mM. The conformations of the FKBP/ligand complex saved along multiple trajectories (50 runs at 310 K for each ligand) are grouped according to a set of intermolecular distances into nodes of a network, and the direct transitions between them are the links. The network analysis reveals that the bound state consists of several subbasins, i.e., binding modes characterized by distinct intermolecular hydrogen bonds and hydrophobic contacts. The dissociation kinetics show a simple (i.e., single-exponential) time dependence because the unbinding barrier is much higher than the barriers between subbasins in the bound state. The unbinding transition state is made up of heterogeneous positions and orientations of the ligand in the FKBP active site, which correspond to multiple pathways of dissociation. For the six small ligands of FKBP, the weaker the binding affinity the closer to the bound state (along the intermolecular distance) are the transition state structures, which is a new manifestation of Hammond behavior. Experimental approaches to the study of fragment binding to proteins have limitations in temporal and spatial resolution. Our network analysis of the unbinding simulations of small inhibitors from an enzyme paints a clear picture of the free energy landscape (both thermodynamics and kinetics) of ligand unbinding.     

The FKBP38 catalytic domain binds to Bcl-2 via a charge-sensitive loop

FKBP38 is a regulator of the prosurvival protein Bcl-2, but in the absence of detailed structural insights, the molecular mechanism of the underlying interaction has remained unknown. Here, we report the contact regions between Bcl-2 and the catalytic domain of FKBP38 derived by heteronuclear NMR spectroscopy. The data reveal that a previously identified charge-sensitive loop near the putative active site of FKBP38 is mainly responsible for Bcl-2 binding. The corresponding binding epitope of Bcl-2 could be identified via a peptide library-based membrane assay. Site-directed mutagenesis of the key residues verified the contact sites of this electrostatic protein/protein interaction. The derived structure model of the complex between Bcl-2 and the FKBP38 catalytic domain features both electrostatic and hydrophobic intermolecular contacts and provides a rationale for the regulation of the FKBP38/Bcl-2 interaction by Ca(2+).     

A high-throughput fluorescence chemical denaturation assay as a general screen for protein–ligand binding

Chemical denaturation of ligand–protein complexes can provide the basis of a label-free binding assay. Here, we show how the technique can be used as a sensitive/affordable screen of potential ligands from a pool of lead drug variants. To demonstrate, we characterized the binding of polyketide ligands based on the mTOR inhibitor rapamycin to the cellular immunophilin FKBP12. This used the intrinsic fluorescence of the protein to monitor the chemical denaturation of each FKBP12–ligand complex. The assay was then successfully modified to a 96-well plate-based screen. Both formats were able to differentiate binding affinities across a wide dynamic range.     

An automated method for dynamic ligand design

An automated method for the dynamic ligand design (DLD) for a binding site of known structure is described. The method can be used for the creation of de novo ligands and for the modification of existing ligands. The binding site is saturated with atoms (sp³ carbon atoms in the present implementation) that form molecules under the influence of a potential function that joins atoms to each other with the correct stereochemistry. The resulting molecules are linked to precomputed functional group minimum energy positions in the binding site. The generalized potential function allows atoms to sample a continuous parameter space that includes the Cartesian coordinates and their occupancy and type, e.g., the method allows change of an sp³ carbon into an sp² carbon or oxygen. A parameter space formulated in this way can then be sampled and optimized by a variety of methods. In this work, molecules are generated by use of a Monte Carlo simulated annealing algorithm. The DLD method is illustrated by its application to the binding site of FK506 binding protein (FKBP), an immunophilin. De novo ligands are designed and modification of the immunosuppressant drug FK506 are suggested. The results demonstrate that the dynamic ligand design approach can automatically construct ligands which complement both the shape and charge distribution of the binding site. © 1995 Wiley-Liss, Inc.     

Switches of hydrogen bonds during ligand–protein association processes determine binding kinetics

Revealing the processes of ligand–protein associations deepens our understanding of molecular recognition and binding kinetics. Hydrogen bonds (H‐bonds) play a crucial role in optimizing ligand–protein interactions and ligand specificity. In addition to the formation of stable H‐bonds in the final bound state, the formation of transient H‐bonds during binding processes contributes binding kinetics that define a ligand as a fast or slow binder, which also affects drug action. However, the effect of forming the transient H‐bonds on the kinetic properties is little understood. Guided by results from coarse‐grained Brownian dynamics simulations, we used classical molecular dynamics simulations in an implicit solvent model and accelerated molecular dynamics simulations in explicit waters to show that the position and distribution of the H‐bond donor or acceptor of a drug result in switching intermolecular and intramolecular H‐bond pairs during ligand recognition processes. We studied two major types of HIV‐1 protease ligands: a fast binder, xk263, and a slow binder, ritonavir. The slow association rate in ritonavir can be attributed to increased flexibility of ritonavir, which yields multistep transitions and stepwise entering patterns and the formation and breaking of complex H‐bond pairs during the binding process. This model suggests the importance of conversions of spatiotemporal H‐bonds during the association of ligands and proteins, which helps in designing inhibitors with preferred binding kinetics. Copyright © 2014 John Wiley & Sons, Ltd.     

Absolute binding free energy calculations using molecular dynamics simulations with restraining potentials

The absolute (standard) binding free energy of eight FK506-related ligands to FKBP12 is calculated using free energy perturbation molecular dynamics (FEP/MD) simulations with explicit solvent. A number of features are implemented to improve the accuracy and enhance the convergence of the calculations. First, the absolute binding free energy is decomposed into sequential steps during which the ligand-surrounding interactions as well as various biasing potentials restraining the translation, orientation, and conformation of the ligand are turned "on" and "off." Second, sampling of the ligand conformation is enforced by a restraining potential based on the root mean-square deviation relative to the bound state conformation. The effect of all the restraining potentials is rigorously unbiased, and it is shown explicitly that the final results are independent of all artificial restraints. Third, the repulsive and dispersive free energy contribution arising from the Lennard-Jones interactions of the ligand with its surrounding (protein and solvent) is calculated using the Weeks-Chandler-Andersen separation. This separation also improves convergence of the FEP/MD calculations. Fourth, to decrease the computational cost, only a small number of atoms in the vicinity of the binding site are simulated explicitly, while all the influence of the remaining atoms is incorporated implicitly using the generalized solvent boundary potential (GSBP) method. With GSBP, the size of the simulated FKBP12/ligand systems is significantly reduced, from approximately 25,000 to 2500. The computations are very efficient and the statistical error is small ( approximately 1 kcal/mol). The calculated binding free energies are generally in good agreement with available experimental data and previous calculations (within approximately 2 kcal/mol). The present results indicate that a strategy based on FEP/MD simulations of a reduced GSBP atomic model sampled with conformational, translational, and orientational restraining potentials can be computationally inexpensive and accurate.     

Template-based docking of a prolactin receptor proline-rich motif octapeptide to FKBP12: implications for cytokine receptor signaling.

A conserved proline-rich motif (PRM) in the cytoplasmic domain of cytokine receptors has been suggested to be a signaling switch regulated by the action of the FK506 binding protein (FKBP) family of peptidylprolyl isomerases (O'Neal KD, Yu-Lee LY, Shearer WT, 1995, Ann NY Acad Sci 766:282-284). We have docked the prolactin receptor PRM (Ile1-Phe2-Pro3-Pro4-Val5-Pro6-Gly7-Pro8) to the ligand binding site of FKBP12. The procedure involved conformational search restricted by NMR restraints (O'Neal KD et al., 1996, Biochem J 315:833-844), energy minimization of the octapeptide conformers so obtained, template-based docking of a selected conformer to FKBP12, and energy refinement of the resulting complex. The template used was the crystal structure of a cyclic FK506-peptide hybrid bound to FKBP12. Val5-Pro6 of the PRM was taken to be the biologically relevant Xaa-Pro bond. The docked conformer is stabilized by two intramolecular hydrogen bonds, N7H7-->O4 and N2H2-->O8, and two intermolecular ones, Ile56; N-H-->O = C:Pro6 and Tyr82:O-H-->O = C:Gly7. This conformer features a Type I beta-turn and has extensive hydrophobic contacts with the FKBP12 binding surface. The observed interactions support the hypothesis that FKBP12 catalyzes cis-trans isomerization in the PRM when it is part of the longer cytoplasmic domain of a cytokine receptor, and suggest a significant role for the PRM in signal transduction.     

Charged surface residues of FKBP12 participate in formation of the FKBP12-FK506-calcineurin complex.

The mechanism of FK506 immunosuppression has been proposed to proceed by formation of a tight-binding complex with the intracellular 12-kDa FK506-binding protein (FKBP12). The FK506-FKBP12 complex then acts as a specific high-affinity inhibitor of the intracellular protein phosphatase PP2B (calcineurin), interrupting downstream dephosphorylation events required for T-cell activation. Site-directed mutagenesis of many of the surface residues of FKBP12 has no significant effect on its affinity for calcineurin. We have identified, however, three FKBP12 surface residues (Asp-37, Arg-42, and His-87) proximal to a solvent-exposed segment of bound FK506 that may be direct contacts in the calcineurin complex. Site-directed mutagenesis of two of these residues decreases the affinity of FKBP12-FK506 for calcineurin (Ki) from 6 nM for wild-type FKBP12 to 3.7 microM for a R42K/H87V double mutant, without affecting the peptidylprolyl isomerase activity or FK506 affinity of the mutant protein. These FKBP12 mutations along with several substitutions on FK506 known to affect calcineurin binding form a roughly 100-A2 region of the FKBP12-FK506 complex surface that is likely to be within the calcineurin binding site.     

Rapid protein-ligand docking using soft modes from molecular dynamics simulations to account for protein deformability: binding of FK506 to FKBP

Most current docking methods to identify possible ligands and putative binding sites on a receptor molecule assume a rigid receptor structure to allow virtual screening of large ligand databases. However, binding of a ligand can lead to changes in the receptor protein conformation that are sterically necessary to accommodate a bound ligand. An approach is presented that allows relaxation of the protein conformation in precalculated soft flexible degrees of freedom during ligand-receptor docking. For the immunosuppressant FK506-binding protein FKBP, the soft flexible modes are extracted as principal components of motion from a molecular dynamics simulation. A simple penalty function for deformations in the soft flexible mode is used to limit receptor protein deformations during docking that avoids a costly recalculation of the receptor energy by summing over all receptor atom pairs at each step. Rigid docking of the FK506 ligand binding to an unbound FKBP conformation failed to identify a geometry close to experiment as favorable binding site. In contrast, inclusion of the flexible soft modes during systematic docking runs selected a binding geometry close to experiment as lowest energy conformation. This has been achieved at a modest increase of computational cost compared to rigid docking. The approach could provide a computationally efficient way to approximately account for receptor flexibility during docking of large numbers of putative ligands and putative docking geometries.     

X-ray structures of small ligand-FKBP complexes provide an estimate for hydrophobic interaction energies

A new crystal form of native FK506 binding protein (FKBP) has been obtained which has proved useful in ligand binding studies. Three different small molecule ligand complexes and the native enzyme have been determined at higher resolution than 2.0 A. Dissociation constants of the related small molecule ligands vary from 20 mM for dimethylsulphoxide to 200 microM for tetrahydrothiophene 1-oxide. Comparison of the four available crystal structures shows that the protein structures are identical to within experimental error, but there are differences in the water structure in the active site. Analysis of the calculated buried surface areas of these related ligands provides an estimated van der Waals contribution to the binding energy of -0.5 kJ/A(2) for non-polar interactions between ligand and protein.     

Molecular docking using surface complementarity

A method is described to dock a ligand into a binding site in a protein on the basis of the complementarity of the inter-molecular atomic contacts. Docking is performed by maximization of a complementarity function that is dependent on atomic contact surface area and the chemical properties of the contacting atoms. The generality and simplicity of the complementarity function ensure that a wide range of chemical structures can be handled. The ligand and the protein are treated as rigid bodies, but displacement of a small number of residues lining the ligand binding site can be taken into account. The method can assist in the design of improved ligands by indicating what changes in complementarity may occur as a result of the substitution of an atom in the ligand. The capabilities of the method are demonstrated by application to 14 protein–ligand complexes of known crystal structure. © 1996 Wiley Liss, Inc.     

The Utility of FK506-Binding Protein as a Fusion Partner in Scintillation Proximity Assays: Application to SH2 Domains

Methodology has been developed which gives a specific measure of the interaction of an SH2 domain with a phosphopeptide ligand using scintillation proximity assay (SPA) technology. Recombinant SH2 domains were expressed from a T7 RNA polymerase-based vector inEscherichia colias fusions to the C-terminus of the FK506-binding protein (FKBP) and purified from freeze-thaw lysates in high yield by affinity chromatography using immobilized phosphopeptides. For binding assays the phosphopeptide ligands were synthesized with a biotin tag and the FKBP fusion proteins were noncovalently radiolabeled with commercially available [³H]dihydroFK506. Complexes of tritiated SH2 fusion protein and biotinyl-phosphopeptide were then captured on streptavidin-coated SPA beads and counted. The modular protocol is an equilibrium technique that does not employ washing steps or specialized radiochemical syntheses required in other binding assays. The utility of the assay has been demonstrated in an examination of the ligand specificity of the SH2 domains of the tyrosine kinases ZAP70, Syk, and Lck. The methodology is potentially generalizable to any receptor–ligand interaction in which one component can be expressed as a fusion partner with FKBP and the other component can be captured on a SPA bead.     

A yeast sensor of ligand binding.

We describe a biosensor that reports the binding of small-molecule ligands to proteins as changes in growth of temperature-sensitive yeast. The yeast strains lack dihydrofolate reductase (DHFR) and are complemented by mouse DHFR containing a ligand-binding domain inserted in a flexible loop. Yeast strains expressing two ligand-binding domain fusions, FKBP12-DHFR and estrogen receptor-alpha (ERalpha)-DHFR, show increased growth in the presence of their corresponding ligands. We used this sensor to identify mutations in residues of ERalpha important for ligand binding, as well as mutations generally affecting protein activity or expression. We also tested the sensor against a chemical array to identify ligands that bind to FKBP12 or ERalpha. The ERalpha sensor was able to discriminate among estrogen analogs, showing different degrees of growth for the analogs that correlated with their relative binding affinities (RBAs). This growth assay provides a simple and inexpensive method to select novel ligands and ligand-binding domains.     

An autonomously folding beta-hairpin derived from the human YAP65 WW domain: attempts to define a minimum ligand-binding motif

WW domains are broadly distributed among natural proteins; these modules play a role in bringing specific proteins together. The ligands recognized by WW domains are short segments rich in proline residues. We have tried to identify the minimum substructure within a WW domain that is required for ligand binding. WW domains typically comprise ca. 40 residues and fold to a three-stranded beta-sheet. Structural data for several WW domain/ligand complexes suggest that most or all of the intermolecular contacts involve beta-strands 2 and 3. We have developed a 16-residue peptide that folds to a beta-hairpin conformation that appears to mimic beta-strands 2 and 3 of the human YAP65 WW domain, but this peptide does not bind to known ligands. Thus, the minimum binding domain is larger than the latter two strands of the WW domain beta-sheet.     

Regulation of AMPA receptor gating by ligand binding core dimers

Ionotropic glutamate receptors are tetramers, the isolated ligand binding cores of which assemble as dimers. Previous work on nondesensitizing AMPA receptor mutants, which combined crystallography, ultracentrifugation, and patch-clamp recording, showed that dimer formation by the ligand binding cores is required for activation of ion channel gating by agonists. To define the mechanisms responsible for stabilization of dimer assembly in native AMPA receptors, contacts between the adjacent ligand binding cores were individually targeted by amino acid substitutions, using the GluR2 crystal structure as a guide to design mutants. We show that disruption of a salt bridge, hydrogen bond network, and intermolecular van der Waals contacts between helices D and J in adjacent ligand binding cores greatly accelerates desensitization. Conservation of these contacts in AMPA and kainate receptors indicates that they are important determinants of dimer stability and that the dimer interface is a key structural element in the gating mechanism of these glutamate receptor families.     

Current noise of a protein-selective biological nanopore

1/f current noise is ubiquitous in protein pores, porins, and channels. We have previously shown that a protein-selective biological nanopore with an external protein receptor can function as a 1/f noise generator when a high-affinity protein ligand is reversibly captured by the receptor. Here, we demonstrate that the binding affinity and concentration of the ligand are key determinants for the nature of current noise. For example, 1/f was absent when a protein ligand was reversibly captured at a much lower concentration than its equilibrium dissociation constant against the receptor. Furthermore, we also analyzed the composite current noise that resulted from mixtures of low-affinity and high-affinity ligands against the same receptor. This study highlights the significance of protein recognition events in the current noise fluctuations across biological membranes.     

Validation of an Allosteric Binding Site of Src Kinase Identified by Unbiased Ligand Binding Simulations

Allostery plays a primary role in regulating protein activity, making it an important mechanism in human disease and drug discovery. Identifying allosteric regulatory sites to explore their biological significance and therapeutic potential is invaluable to drug discovery; however, identification remains a challenge. Allosteric sites are often “cryptic” without clear geometric or chemical features. Since allosteric regulatory sites are often less conserved in protein kinases than the orthosteric ATP binding site, allosteric ligands are commonly more specific than ATP competitive inhibitors. We present a generalizable computational protocol to predict allosteric ligand binding sites based on unbiased ligand binding simulation trajectories. We demonstrate the feasibility of this protocol by revisiting our previously published ligand binding simulations using the first identified viral proto-oncogene, Src kinase, as a model system. The binding paths for kinase inhibitor PP1 uncovered three metastable intermediate states before binding the high-affinity ATP-binding pocket, revealing two previously known allosteric sites and one novel site. Herein, we validate the novel site using a combination of virtual screening and experimental assays to identify a V-type allosteric small-molecule inhibitor that targets this novel site with specificity for Src over closely related kinases. This study provides a proof-of-concept for employing unbiased ligand binding simulations to identify cryptic allosteric binding sites and is widely applicable to other protein–ligand systems.     

Probing molecular structure of dioxygen reduction site of bacterial quinol oxidases through ligand binding to the redox metal centers

Cytochromes bo and bd are structurally unrelated terminal ubiquinol oxidases in the aerobic respiratory chain of Escherichia coli. The high-spin heme o-CuB binuclear center serves as the dioxygen reduction site for cytochrome bo, and the heme b595-heme d binuclear center for cytochrome bd. CuB coordinates three histidine ligands and serves as a transient ligand binding site en route to high-spin heme o one-electron donor to the oxy intermediate, and a binding site for bridging ligands like cyanide. In addition, it can protect the dioxygen reduction site through binding of a peroxide ion in the resting state, and connects directly or indirectly Tyr288 and Glu286 to carry out redox-driven proton pumping in the catalytic cycle. Contrary, heme b595 of cytochrome bd participate a similar role to CuB in ligand binding and dioxygen reduction but cannot perform such versatile roles because of its rigid structure.     

Transient aggregation and stable dimerization induced by introducing an Alzheimer sequence into a water-soluble protein

Transient contacts between denatured polypeptide chains are likely to play an important part in the initial stages of protein aggregation and fibrillation. To analyze the nature of such contacts, we have carried out a protein engineering study of the 102-residue protein U1A, which aggregates transiently in the wild-type form during refolding from the guanidinium chloride-denatured state. We have prepared a series of mutants with increased aggregation tendencies by increasing the homology between two beta-strands of U1A and the Alzheimer peptide (beta-AP). These mutants undergo transient aggregation during refolding, as measured by concentration dependence, double-jump experiments, and binding of ANS, a probe for exposed hydrophobic patches on protein surfaces. The propensity to aggregate increases with increasing homology to beta-AP. Further, the degree of transient ANS binding correlates reasonably well with the structural parameters recently shown to play a role in the fibrillation of natively unfolded proteins. Two mutants highly prone to transient aggregation, U1A-J and U1A-G, were also studied by NMR. Secondary structural elements of the U1A-J construct (with lower beta-AP homology) are very similar to those observed in U1A-wt. In contrast, the high-homology construct U1A-G exhibits local unfolding of the C-terminal helix, which packs against the beta-sheet in the wild-type protein. U1A-G is mainly dimeric according to (15)N spin relaxation data, and the dimer interface most likely involves the beta-sheet. Our data suggest that the transient aggregate relies on specific intermolecular interactions mediated by structurally flexible regions and that contacts may be formed in different beta-strand registers.     

154 Protein folding and binding funnels: a common driving force and a common mechanism

Under the free energy landscape theory, both the protein-folding and protein–ligand binding processes are driven by the decrease in total Gibbs free energy of the protein-solvent or protein–ligand-solvent system, which involves the non-complementary changes between the entropy and enthalpy, ultimately leading to a global free energy minimization of these thermodynamic systems (Ji & Liu, 2011; Liu et al., 2012; Yang, Ji & Liu, 2012). In the case of protein folding, the lowering of the system free energy coupled with the gradual reduction in conformational degree of freedom of the folding intermediates determines that the shape of the free energy landscape for protein folding must be funnel-like (Dill & Chan, 1997), rather than non-funneled shapes (Ben-Naim, 2012). In the funnel-like free energy landscape, protein folding can be viewed as going down the hill via multiple parallel routes from a vast majority of individual non-native states on surface outside the funnel to the native states located around the bottom of the funnel. The first stage of folding, i.e. the rapid hydrophobic collapse process, is driven by the solvent entropy maximization. Concretely, the water molecules squeeze and sequestrate the hydrophobic amino acid side chains within the interior of the folding intermediates while exposing the polar and electrostatically charged side chains on the intermediate surface so as to minimize the solvent-accessible surface area of the solute and thus, the minimal contacts between the folding intermediates and the water molecules. This will maximize the entropy of the solvent, thus contributing substantially to lowering of the system free energy due to an absolute advantage of the solvent in both quantity and mass (Yang, Ji & Liu, 2012). The resulting molten globule states (Ohgushi & Wada, 1983), within which a few transient secondary structural components and tertiary contacts have been formed but many native contacts or close residue–residue interactions has yet to form, need to be further sculptured into the native states. This is a relatively slow “bottleneck” process because the competitive interactions between protein residues within the folding intermediates and between residues and water molecules may repeat many rounds to accumulate a large enough number of stable noncovalent bonds capable of counteracting the conformational entropy loss of the intermediates, thus putting this bottleneck stage under the enthalpy control (i.e. negative enthalpy change), contributing further to the lowering of the system free energy. Although the protein–ligand association occurs around the rugged bottom of the free energy landscape, the exclusion of water from the binding interfaces and the formation of noncovalent bonds between the two partners can still lower the system free energy. In conjunction with the loss of the rotational and translational degrees of freedom of the two partners as well as the loss of the conformational entropy of the protein, these processes could merge, downwards expand, and further narrow the free energy wells within which the protein–ligand binding process takes place, thereby making them look like a funnel, which we term the binding funnel. In this funnel, the free energy downhill process follows a similar paradigm to the protein-folding process. For example, if the initial collisions/contacts occur between the properly complementary interfaces of the protein and ligand, a large amount of water molecules (which usually form a water network around the solute surface) will be displaced to suit the need for maximizing the solvent entropy. This process is similar to that of the hydrophobic collapse during protein folding, resulting in a loosely associated protein–ligand complex that needs also to be further adapted into a tight complex, i.e. the second step which is mainly driven by the negative enthalpy change through intermolecular competitive interactions to gradually accumulate the noncovalent bonds and ultimately, to stabilize the complex at a tightly bound state. Taken together, we conclude that whether in the protein-folding or in the protein–ligand binding process, both the entropy-driven first step and the enthalpy-driven second step contribute to the lowering of the system free energy, resulting in the funnel-like folding or binding free energy landscape.