Adrenomedullin 2/intermedin is a slow off-rate, long-acting endogenous agonist of the adrenomedullin2 G protein-coupled receptor.

Adrenomedullin 2/intermedin (AM2/IMD), adrenomedullin (AM), and calcitonin gene-related peptide (CGRP) have functions in the cardiovascular, lymphatic, and nervous systems by activating three heterodimeric receptors comprising the class B GPCR CLR and a RAMP1, -2, or -3 modulatory subunit. CGRP and AM prefer the RAMP1 and RAMP2/3 complexes, respectively, whereas AM2/IMD is thought to be relatively nonselective. Accordingly, AM2/IMD exhibits overlapping actions with CGRP and AM, so the rationale for this third agonist for the CLR-RAMP complexes is unclear. Here, we report that AM2/IMD is kinetically selective for CLR-RAMP3, known as the AM2R, and we define the structural basis for its distinct kinetics. In live cell biosensor assays, AM2/IMD-AM2R elicited longer-duration cAMP signaling than the other peptide-receptor combinations. AM2/IMD and AM bound the AM2R with similar equilibrium affinities, but AM2/IMD had a slower off-rate and longer receptor residence time, thus explaining its prolonged signaling capacity. Peptide and receptor chimeras and mutagenesis were used to map the regions responsible for the distinct binding and signaling kinetics to the AM2/IMD mid-region and the RAMP3 extracellular domain (ECD). Molecular dynamics simulations revealed how the former forms stable interactions at the CLR ECD-transmembrane domain interface and how the latter augments the CLR ECD binding pocket to anchor the AM2/IMD C terminus. These strong binding components only combine in the AM2R. Our findings uncover AM2/IMD-AM2R as a cognate pair with unique temporal features, reveal how AM2/IMD and RAMP3 collaborate to shape CLR signaling, and have significant implications for AM2/IMD biology.

The Small Molecule Antagonist KCI807 Disrupts Association of the Amino-Terminal Domain of the Androgen Receptor with ELK1 by Modulating the Adjacent DNA Binding Domain.

The androgen receptor (AR) is a crucial coactivator of ELK1 for prostate cancer (PCa) growth, associating with ELK1 through two peptide segments (358-457 and 514-557) within the amino-terminal domain (NTD) of AR. The small-molecule antagonist 5-hydroxy-2-(3-hydroxyphenyl)chromen-4-one (KCI807) binds to AR, blocking ELK1 binding and inhibiting PCa growth. We investigated the mode of interaction of KCI807 with AR using systematic mutagenesis coupled with ELK1 coactivation assays, testing polypeptide binding and Raman spectroscopy. In full-length AR, deletion of neither ELK1 binding segment affected sensitivity of residual ELK1 coactivation to KCI807. Although the NTD is sufficient for association of AR with ELK1, interaction of the isolated NTD with ELK1 was insensitive to KCI807. In contrast, coactivation of ELK1 by the AR-V7 splice variant, comprising the NTD and the DNA binding domain (DBD), was sensitive to KCI807. Deletions and point mutations within DBD segment 558-595, adjacent to the NTD, interfered with coactivation of ELK1, and residual ELK1 coactivation by the mutants was insensitive to KCI807. In a glutathione S-transferase pull-down assay, KCI807 inhibited ELK1 binding to an AR polypeptide that included the two ELK1 binding segments and the DBD but did not affect ELK1 binding to a similar AR segment that lacked the sequence downstream of residue 566. Raman spectroscopy detected KCI807-induced conformational change in the DBD. The data point to a putative KCI807 binding pocket within the crystal structure of the DBD and indicate that either mutations or binding of KCI807 at this site will induce conformational changes that disrupt ELK1 binding to the NTD. SIGNIFICANCE STATEMENT: The small-molecule antagonist KCI807 disrupts association of the androgen receptor (AR) with ELK1, serving as a prototype for the development of small molecules for a novel type of therapeutic intervention in drug-resistant prostate cancer. This study provides basic information needed for rational KCI807-based drug design by identifying a putative binding pocket in the DNA binding domain of AR through which KCI807 modulates the amino-terminal domain to inhibit ELK1 binding.

A biosensor of protein foldedness identifies increased "holdase" activity of chaperones in the nucleus following increased cytosolic protein aggregation.

Chaperones and other quality control machinery guard proteins from inappropriate aggregation, which is a hallmark of neurodegenerative diseases. However, how the systems that regulate the "foldedness" of the proteome remain buffered under stress conditions and in different cellular compartments remains incompletely understood. In this study, we applied a FRET-based strategy to explore how well quality control machinery protects against the misfolding and aggregation of "bait" biosensor proteins, made from the prokaryotic ribonuclease barnase, in the nucleus and cytosol of human embryonic kidney 293T cells. We found that those barnase biosensors were prone to misfolding, were less engaged by quality control machinery, and more prone to inappropriate aggregation in the nucleus as compared with the cytosol, and that these effects could be regulated by chaperone Hsp70-related machinery. Furthermore, aggregation of mutant huntingtin exon 1 protein (Httex1) in the cytosol appeared to outcompete and thus prevented the engagement of quality control machinery with the biosensor in the cytosol. This effect correlated with reduced levels of DNAJB1 and HSPA1A chaperones in the cell outside those sequestered to the aggregates, particularly in the nucleus. Unexpectedly, we found Httex1 aggregation also increased the apparent engagement of the barnase biosensor with quality control machinery in the nucleus suggesting an independent implementation of "holdase" activity of chaperones other than DNAJB1 and HSPA1A. Collectively, these results suggest that proteostasis stress can trigger a rebalancing of chaperone abundance in different subcellular compartments through a dynamic network involving different chaperone-client interactions.

How Robust Is the Ligand Binding Transition State?

For many drug targets, it has been shown that the kinetics of drug binding (e.g., on rate and off rate) is more predictive of drug efficacy than thermodynamic quantities alone. This motivates the development of predictive computational models that can be used to optimize compounds on the basis of their kinetics. The structural details underpinning these computational models are found not only in the bound state but also in the short-lived ligand binding transition states. Although transition states cannot be directly observed experimentally due to their extremely short lifetimes, recent successes have demonstrated that modeling the ligand binding transition state is possible with the help of enhanced sampling molecular dynamics methods. Previously, we generated unbinding paths for an inhibitor of soluble epoxide hydrolase (sEH) with a residence time of 11 min. Here, we computationally modeled unbinding events with the weighted ensemble method REVO (resampling of ensembles by variation optimization) for five additional inhibitors of sEH with residence times ranging from 14.25 to 31.75 min, with average prediction accuracy within an order of magnitude. The unbinding ensembles are analyzed in detail, focusing on features of the ligand binding transition state ensembles (TSEs). We find that ligands with similar bound poses can show significant differences in their ligand binding TSEs, in terms of their spatial distribution and protein-ligand interactions. However, we also find similarities across the TSEs when examining more general features such as ligand degrees of freedom. Together these findings show significant challenges for rational, kinetics-based drug design.

Quality over quantity: Sampling high probability rare events with the weighted ensemble algorithm.

The prediction of (un)binding rates and free energies is of great significance to the drug design process. Although many enhanced sampling algorithms and approaches have been developed, there is not yet a reliable workflow to predict these quantities. Previously we have shown that free energies and transition rates can be calculated by directly simulating the binding and unbinding processes with our variant of the WE algorithm "Resampling of Ensembles by Variation Optimization", or "REVO". Here, we calculate binding free energies retrospectively for three SAMPL6 host-guest systems and prospectively for a SAMPL9 system to test a modification of REVO that restricts its cloning behavior in quasi-unbound states. Specifically, trajectories cannot clone if they meet a physical requirement that represents a high likelihood of unbinding, which in the case of this work is a center-of-mass to center-of-mass distance. The overall effect of this change was difficult to predict, as it results in fewer unbinding events each of which with a much higher statistical weight. For all four systems tested, this new strategy produced either more accurate unbinding free energies or more consistent results between simulations than the standard REVO algorithm. This approach is highly flexible, and any feature of interest for a system can be used to determine cloning eligibility. These findings thus constitute an important improvement in the calculation of transition rates and binding free energies with the weighted ensemble method.

Membrane-Mediated Ligand Unbinding of the PK-11195 Ligand from TSPO.

The translocator protein (TSPO), previously known as the peripheral benzodiazepine receptor, is of longstanding medical interest as both a biomarker for neuroinjury and a potential drug target for neuroinflammation and other disorders. Recently, it was shown that ligand residence time is a key factor determining steroidogenic efficacy of TSPO-binding compounds. This spurs interest in simulations of (un)binding pathways of TSPO ligands, which could reveal the molecular interactions governing ligand residence time. In this study, we use a weighted ensemble algorithm to determine the unbinding pathway for different poses of PK-11195, a TSPO ligand used in neuroimaging. In contrast with previous studies, our results show that PK-11195 does not dissociate directly into the solvent but instead dissociates via the lipid membrane by going between the transmembrane helices. We analyze this path ensemble in detail, constructing descriptors that can facilitate a general understanding of membrane-mediated ligand binding. We construct a set of Markov state models augmented with additional straightforward simulations to determine pose-specific ligand residence times. Together, we combine over 40 µs of trajectory data to form a coherent picture of the ligand binding landscape. We find that multiple starting poses yield residence times that roughly agree with the experimental quantity. The ligand binding transition states predicted by these Markov state models occur when PK-11195 is already in the membrane and involves only minimal ligand-protein interactions. This has implications for the design of new long-residence-time TSPO ligands.

ClassicalGSG: Prediction of log P using classical molecular force fields and geometric scattering for graphs.

This work examines methods for predicting the partition coefficient (log P) for a dataset of small molecules. Here, we use atomic attributes such as radius and partial charge, which are typically used as force field parameters in classical molecular dynamics simulations. These atomic attributes are transformed into index-invariant molecular features using a recently developed method called geometric scattering for graphs (GSG). We call this approach "ClassicalGSG" and examine its performance under a broad range of conditions and hyperparameters. We train ClassicalGSG log P predictors with neural networks using 10,722 molecules from the OpenChem dataset and apply them to predict the log P values from four independent test sets. The ClassicalGSG method's performance is compared to a baseline model that employs graph convolutional networks. Our results show that the best prediction accuracies are obtained using atomic attributes generated with the CHARMM generalized force field and 2D molecular structures.

Identification and Characterization of Two Structurally Related Dipeptides that Enhance Catalytic Efficiency of Neurolysin.

Neurolysin (Nln) is a recently recognized endogenous mechanism functioning to preserve the brain from ischemic injury. To further understand the pathophysiological function of this peptidase in stroke and other neurologic disorders, the present study was designed to identify small molecule activators of Nln. Using a computational approach, the structure of Nln was explored, which was followed by docking and in silico screening of ∼140,000 molecules from the National Cancer Institute Developmental Therapeutics Program database. Top ranking compounds were evaluated in an Nln enzymatic assay, and two hit histidine-dipeptides were further studied in detail. The identified dipeptides enhanced the rate of synthetic substrate hydrolysis by recombinant (human and rat) and mouse brain-purified Nln in a concentration-dependent manner (micromolar A50 and Amax ≥ 300%) but had negligible effect on activity of closely related peptidases. Both dipeptides also enhanced hydrolysis of Nln endogenous substrates neurotensin, angiotensin I, and bradykinin and increased efficiency of the synthetic substrate hydrolysis (Vmax/Km ratio) in a concentration-dependent manner. The dipeptides and competitive inhibitor dynorphin A (1-13) did not affect each other's affinity for Nln, suggesting differing nature of their respective binding sites. Lastly, drug affinity responsive target stability (DARTS) and differential scanning fluorimetry (DSF) assays confirmed concentration-dependent interaction of Nln with the activator molecule. This is the first study demonstrating that Nln activity can be enhanced by small molecules, although the peptidic nature and low potency of the activators limit their application. The identified dipeptides provide a chemical scaffold to develop high-potency, drug-like molecules as research tools and potential drug leads. SIGNIFICANCE STATEMENT: This study describes discovery of two molecules that selectively enhance activity of peptidase Nln-a newly recognized cerebroprotective mechanism in the poststroke brain. The identified molecules will serve as a chemical scaffold for development of drug-like molecules to further study Nln and may become lead structures for a new class of drugs. In addition, our conceptual and methodological framework and research findings might be used for other peptidases and enzymes, the activation of which bears therapeutic potential.

Predicting partition coefficients for the SAMPL7 physical property challenge using the ClassicalGSG method.

The prediction of [Formula: see text] values is one part of the statistical assessment of the modeling of proteins and ligands (SAMPL) blind challenges. Here, we use a molecular graph representation method called Geometric Scattering for Graphs (GSG) to transform atomic attributes to molecular features. The atomic attributes used here are parameters from classical molecular force fields including partial charges and Lennard-Jones interaction parameters. The molecular features from GSG are used as inputs to neural networks that are trained using a "master" dataset comprised of over 41,000 unique [Formula: see text] values. The specific molecular targets in the SAMPL7 [Formula: see text] prediction challenge were unique in that they all contained a sulfonyl moeity. This motivated a set of ClassicalGSG submissions where predictors were trained on different subsets of the master dataset that are filtered according to chemical types and/or the presence of the sulfonyl moeity. We find that our ranked prediction obtained 5th place with an RMSE of 0.77 [Formula: see text] units and an MAE of 0.62, while one of our non-ranked predictions achieved first place among all submissions with an RMSE of 0.55 and an MAE of 0.44. After the conclusion of the challenge we also examined the performance of open-source force field parameters that allow for an end-to-end [Formula: see text] predictor model: General AMBER Force Field (GAFF), Universal Force Field (UFF), Merck Molecular Force Field 94 (MMFF94) and Ghemical. We find that ClassicalGSG models trained with atomic attributes from MMFF94 can yield more accurate predictions compared to those trained with CGenFF atomic attributes.

The SAMPL6 SAMPLing challenge: assessing the reliability and efficiency of binding free energy calculations.

Approaches for computing small molecule binding free energies based on molecular simulations are now regularly being employed by academic and industry practitioners to study receptor-ligand systems and prioritize the synthesis of small molecules for ligand design. Given the variety of methods and implementations available, it is natural to ask how the convergence rates and final predictions of these methods compare. In this study, we describe the concept and results for the SAMPL6 SAMPLing challenge, the first challenge from the SAMPL series focusing on the assessment of convergence properties and reproducibility of binding free energy methodologies. We provided parameter files, partial charges, and multiple initial geometries for two octa-acid (OA) and one cucurbit[8]uril (CB8) host-guest systems. Participants submitted binding free energy predictions as a function of the number of force and energy evaluations for seven different alchemical and physical-pathway (i.e., potential of mean force and weighted ensemble of trajectories) methodologies implemented with the GROMACS, AMBER, NAMD, or OpenMM simulation engines. To rank the methods, we developed an efficiency statistic based on bias and variance of the free energy estimates. For the two small OA binders, the free energy estimates computed with alchemical and potential of mean force approaches show relatively similar variance and bias as a function of the number of energy/force evaluations, with the attach-pull-release (APR), GROMACS expanded ensemble, and NAMD double decoupling submissions obtaining the greatest efficiency. The differences between the methods increase when analyzing the CB8-quinine system, where both the guest size and correlation times for system dynamics are greater. For this system, nonequilibrium switching (GROMACS/NS-DS/SB) obtained the overall highest efficiency. Surprisingly, the results suggest that specifying force field parameters and partial charges is insufficient to generally ensure reproducibility, and we observe differences between seemingly converged predictions ranging approximately from 0.3 to 1.0 kcal/mol, even with almost identical simulations parameters and system setup (e.g., Lennard-Jones cutoff, ionic composition). Further work will be required to completely identify the exact source of these discrepancies. Among the conclusions emerging from the data, we found that Hamiltonian replica exchange-while displaying very small variance-can be affected by a slowly-decaying bias that depends on the initial population of the replicas, that bidirectional estimators are significantly more efficient than unidirectional estimators for nonequilibrium free energy calculations for systems considered, and that the Berendsen barostat introduces non-negligible artifacts in expanded ensemble simulations.

Enhanced Jarzynski free energy calculations using weighted ensemble.

The free energy of transitions between stable states is the key thermodynamic quantity that governs the relative probabilities of the forward and reverse reactions and the ratio of state probabilities at equilibrium. The binding free energy of a drug and its receptor is of particular interest, as it serves as an optimization function for drug design. Over the years, many computational methods have been developed to calculate binding free energies, and while many of these methods have a long history, issues such as convergence of free energy estimates and the projection of a binding process onto order parameters remain. Over 20 years ago, the Jarzynski equality was derived with the promise to calculate equilibrium free energies by measuring the work applied to short nonequilibrium trajectories. However, these calculations were found to be dominated by trajectories with low applied work that occur with extremely low probability. Here, we examine the combination of weighted ensemble algorithms with the Jarzynski equality. In this combined method, an ensemble of nonequilibrium trajectories are run in parallel, and cloning and merging operations are used to preferentially sample low-work trajectories that dominate the free energy calculations. Two additional methods are also examined: (i) a novel weighted ensemble resampler that samples trajectories directly according to their importance to the work of work and (ii) the diffusion Monte Carlo method using the applied work as the selection potential. We thoroughly examine both the accuracy and efficiency of unbinding free energy calculations for a series of model Lennard-Jones atom pairs with interaction strengths ranging from 2 kcal/mol to 20 kcal/mol. We find that weighted ensemble calculations can more efficiently determine accurate binding free energies, especially for deeper Lennard-Jones well depths.

Selectivity, ligand deconstruction, and cellular activity analysis of a BPTF bromodomain inhibitor.

Bromodomain and PHD finger containing protein transcription factor (BPTF) is an epigenetic protein involved in chromatin remodelling and is a potential anticancer target. The BPTF bromodomain has one reported small molecule inhibitor (AU1, rac-1). Here, advances made on the structure-activity relationship of a BPTF bromodomain ligand are reported using a combination of experimental and molecular dynamics simulations leading to the active enatiomer (S)-1. Additionally, a ligand deconstruction analysis was conducted to characterize important pharmacophores for engaging the BPTF bromodomain. These studies have been enabled by a protein-based fluorine NMR approach, highlighting the versatility of the method for selectivity, ligand deconstruction, and ligand binding. To enable future analysis of biological activity, cell growth analyses in a panel of cancer cell lines were carried out using CRISPR-Cas9 and (S)-1 to identify cell-based model systems that are sensitive to BPTF inhibition.

REVO: Resampling of ensembles by variation optimization.

Conventional molecular dynamics simulations are incapable of sampling many important interactions in biomolecular systems due to their high dimensionality and rough energy landscapes. To observe rare events and calculate transition rates in these systems, enhanced sampling is a necessity. In particular, the study of ligand-protein interactions necessitates a diverse ensemble of protein conformations and transition states, and for many systems, this occurs on prohibitively long time scales. Previous strategies such as WExplore that can be used to determine these types of ensembles are hindered by problems related to the regioning of conformational space. Here, we propose a novel, regionless, enhanced sampling method that is based on the weighted ensemble framework. In this method, a value referred to as "trajectory variation" is optimized after each cycle through cloning and merging operations. This method allows for a more consistent measurement of observables and broader sampling resulting in the efficient exploration of previously unexplored conformations. We demonstrate the performance of this algorithm with the N-dimensional random walk and the unbinding of the trypsin-benzamidine system. The system is analyzed using conformation space networks, the residence time of benzamidine is confirmed, and a new unbinding pathway for the trypsin-benzamidine system is found. We expect that resampling of ensembles by variation optimization will be a useful general tool to broadly explore free energy landscapes.

Structural insights into lethal contractural syndrome type 3 (LCCS3) caused by a missense mutation of PIP5Kγ.

Signaling molecule phosphatidylinositol 4,5-bisphosphate is produced primarily by phosphatidylinositol 4-phosphate 5-kinase (PIP5K). PIP5K is essential for the development of the human neuronal system, which has been exemplified by a recessive genetic disorder, lethal congenital contractural syndrome type 3, caused by a single aspartate-to-asparagine mutation in the kinase domain of PIP5Kγ. So far, the exact role of this aspartate residue has yet to be elucidated. In this work, we conducted structural, functional and computational studies on a zebrafish PIP5Kα variant with a mutation at the same site. Compared with the structure of the wild-type (WT) protein in the ATP-bound state, the ATP-associating glycine-rich loop of the mutant protein was severely disordered and the temperature factor of ATP was significantly higher. Both observations suggest a greater degree of disorder of the bound ATP, whereas neither the structure of the catalytic site nor the Km toward ATP was substantially affected by the mutation. Microsecond molecular dynamics simulation revealed that negative charge elimination caused by the mutation destabilized the involved hydrogen bonds and affected key electrostatic interactions in the close proximity of ATP. Taken together, our data indicated that the disease-related aspartate residue is a key node in the interaction network crucial for effective ATP binding. This work provides a paradigm of how a subtle but critical structural perturbation caused by a single mutation at the ATP-binding site abolishes the kinase activity, emphasizing that stabilizing substrate in a productive conformational state is crucial for catalysis.

Mapping the Ligand Binding Landscape.

The interaction between a ligand and a protein involves a multitude of conformational states. To achieve a particular deeply bound pose, the ligand must search across a rough free-energy landscape with many metastable minima. Creating maps of the ligand binding landscape is a great challenge, as binding and release events typically occur on timescales that are beyond the reach of molecular simulation. The WExplore enhanced sampling method is well suited to build these maps because it is designed to broadly explore free-energy landscapes and is capable of simulating ligand release pathways that occur on timescales as long as minutes. WExplore also uses only unbiased trajectory segments, allowing for the construction of Markov state models (MSMs) and conformation space networks that combine the results of multiple simulations. Here, we use WExplore to study two bromodomain-inhibitor systems using multiple docked starting poses (Brd4-MS436 and Baz2B-ICR7) and synthesize our results using a series of MSMs using time-lagged independent component analysis. Ranking the starting poses by exit rate agrees with the crystal structure pose in both cases. We also predict the most stable pose using the equilibrium populations from the MSM but find that the prediction is not robust as a function of MSM parameters. The simulated trajectories are synthesized into network models that visualize the entire binding landscape for each system, and we examine transition paths between deeply bound stable states. We find that, on average, transitions between deeply bound states convert through the unbound state 81% of the time, implying a trial-and-error approach to ligand binding. We conclude with a discussion of the implications of this result for both kinetics-based drug discovery and virtual screening pipelines that incorporate molecular dynamics.

Unbiased Molecular Dynamics of 11 min Timescale Drug Unbinding Reveals Transition State Stabilizing Interactions.

Ligand (un)binding kinetics is being recognized as a determinant of drug specificity and efficacy in an increasing number of systems. However, the calculation of kinetics and the simulation of drug unbinding is more difficult than computing thermodynamic quantities, such as binding free energies. Here we present the first full simulations of an unbinding process at pharmacologically relevant time scales (11 min), without the use of biasing forces, detailed prior knowledge, or specialized processors using the weighted ensemble based algorithm, WExplore. These simulations show the inhibitor TPPU unbinding from its enzyme target soluble epoxide hydrolase, which is a clinically relevant target that has attracted interest in kinetics optimization in order to increase efficacy. We make use of conformation space networks that allow us to conceptualize unbinding not just as a linear process, but as a network of interconnected states that connect the bound and unbound states. This allows us to visualize patterns in hydrogen-bonding, solvation, and nonequilibrium free energies, without projection onto progress coordinates. The topology and layout of the network reveal multiple unbinding pathways, and other rare events, such as the reversal of ligand orientation within the binding site. Furthermore, we make a prediction of the transition state ensemble, using transition path theory, and identify protein-ligand interactions which are stabilizing to the transition state. Additionally, we uncover trends in ligand and binding site solvation that corroborate experimental evidence from more classical structure kinetics relationships and generate new questions as to the role of drug modifications in kinetics optimization. Finally, from only 6 µs of simulation time we observed 75 unbinding events from which we calculate a residence time of 42 s, and a standard error range of 23 to 280 s. This nearly encompasses the experimental residence time 11 min (660 s). In addition to the insights to sEH inhibitor unbinding, this study shows that simulations of complex processes on timescales as long as minutes are becoming feasible for more researchers to perform.

Predicting ligand binding affinity using on- and off-rates for the SAMPL6 SAMPLing challenge.

Interest in ligand binding kinetics has been growing rapidly, as it is being discovered in more and more systems that ligand residence time is the crucial factor governing drug efficacy. Many enhanced sampling methods have been developed with the goal of predicting ligand binding rates ([Formula: see text]) and/or ligand unbinding rates ([Formula: see text]) through explicit simulation of ligand binding pathways, and these methods work by very different mechanisms. Although there is not yet a blind challenge for ligand binding kinetics, here we take advantage of experimental measurements and rigorously computed benchmarks to compare estimates of [Formula: see text] calculated as the ratio of two rates: [Formula: see text]. These rates were determined using a new enhanced sampling method based on the weighted ensemble framework that we call "REVO": Reweighting of Ensembles by Variance Optimization. This is a further development of the WExplore enhanced sampling method, in which trajectory cloning and merging steps are guided not by the definition of sampling regions, but by maximizing trajectory variance. Here we obtain estimates of [Formula: see text] and [Formula: see text] that are consistent across multiple simulations, with an average log10-scale standard deviation of 0.28 for on-rates and 0.56 for off-rates, which is well within an order of magnitude and far better than previously observed for previous applications of the WExplore algorithm. Our rank ordering of the three host-guest pairs agrees with the reference calculations, however our predicted [Formula: see text] values were systematically lower than the reference by an average of 4.2 kcal/mol. Using tree network visualizations of the trajectories in the REVO algorithm, and conformation space networks for each system, we analyze the results of our sampling, and hypothesize sources of discrepancy between our [Formula: see text] values and the reference. We also motivate the direct inclusion of [Formula: see text] and [Formula: see text] challenges in future iterations of SAMPL, to further develop the field of ligand binding kinetics prediction and modeling.

Multiple Ligand Unbinding Pathways and Ligand-Induced Destabilization Revealed by WExplore.

We report simulations of full ligand exit pathways for the trypsin-benzamidine system, generated using the sampling technique WExplore. WExplore is able to observe millisecond-scale unbinding events using many nanosecond-scale trajectories that are run without introducing biasing forces. The algorithm generates rare events by dividing the coordinate space into regions, on-the-fly, and balancing computational effort between regions through cloning and merging steps, as in the weighted ensemble method. The averaged exit flux yields a ligand exit rate of 180 µs, which is within an order of magnitude of the experimental value. We obtain broad sampling of ligand exit pathways, and visualize our findings using conformation space networks. The analysis shows three distinct exit channels, two of which are formed through large, rare motions of the loop regions in trypsin. This broad set of ligand-bound poses is then used to investigate general properties of ligand binding: we observe both a direct stabilizing effect of ligand-protein interactions and an indirect destabilizing effect on intraprotein interactions that is induced by the ligand. Significantly, the crystallographic binding poses are distinguished not only because their ligands induce large stabilizing effects, but also because they induce relatively low indirect destabilizations.

Optimal allosteric stabilization sites using contact stabilization analysis.

Proteins can be destabilized by a number of environmental factors such as temperature, pH, and mutation. The ability to subsequently restore function under these conditions by adding small molecule stabilizers, or by introducing disulfide bonds, would be a very powerful tool, but the physical principles that drive this stabilization are not well understood. The first problem lies is in choosing an appropriate binding site or disulfide bond location to best confer stability to the active site and restore function. Here, we present a general framework for predicting which allosteric binding sites correlate with stability in the active site. Using the Karanicolas-Brooks Go-like model, we examine the dynamics of the enzyme ß-glucuronidase using an Umbrella Sampling method to thoroughly sample the conformational landscape. Each intramolecular contact is assigned a score termed a "stabilization factor" that measures its correlation with structural changes in the active site. We have carried out this analysis for three different scaling strengths for the intramolecular contacts, and we examine how the calculated stabilization factors depend on the ensemble of destabilized conformations. We further examine a locally destabilized mutant of ß-glucuronidase that has been characterized experimentally, and show that this brings about local changes in the stabilization factors. We find that the proximity to the active site is not sufficient to determine which contacts can confer active site stability. © 2016 Wiley Periodicals, Inc.

Capturing a Dynamic Chaperone-Substrate Interaction Using NMR-Informed Molecular Modeling.

Chaperones maintain a healthy proteome by preventing aggregation and by aiding in protein folding. Precisely how chaperones influence the conformational properties of their substrates, however, remains unclear. To achieve a detailed description of dynamic chaperone-substrate interactions, we fused site-specific NMR information with coarse-grained simulations. Our model system is the binding and folding of a chaperone substrate, immunity protein 7 (Im7), with the chaperone Spy. We first used an automated procedure in which NMR chemical shifts inform the construction of system-specific force fields that describe each partner individually. The models of the two binding partners are then combined to perform simulations on the chaperone-substrate complex. The binding simulations show excellent agreement with experimental data from multiple biophysical measurements. Upon binding, Im7 interacts with a mixture of hydrophobic and hydrophilic residues on Spy's surface, causing conformational exchange within Im7 to slow down as Im7 folds. Meanwhile, the motion of Spy's flexible loop region increases, allowing for better interaction with different substrate conformations, and helping offset losses in Im7 conformational dynamics that occur upon binding and folding. Spy then preferentially releases Im7 into a well-folded state. Our strategy has enabled a residue-level description of a dynamic chaperone-substrate interaction, improving our understanding of how chaperones facilitate substrate folding. More broadly, we validate our approach using two other binding partners, showing that this approach provides a general platform from which to investigate other flexible biomolecular complexes through the integration of NMR data with efficient computational models.