Potentials of mean force for acetylcholine unbinding from the alpha7 nicotinic acetylcholine receptor ligand-binding domain

Molecular Dynamics Simulation Study on the Structural and Thermodynamic Analysis of Oxidized and Unoxidized Forms of Polyaniline

The conducting polymer polyaniline (PANI) has shown significant interest for the development of electrified membranes (EMs) with superior antifouling characteristics. However, the blending and doping of PANI with other polymers and nanomaterials highly influence the properties of the membrane surface. PANI exists in two forms: oxidized, known as emeraldine salt (ES), and unoxidized, referred to as emeraldine base (EB). Therefore, understanding the different forms of PANI and the variations between the oxidized and unoxidized forms along the length of the polymer chain is intriguing. In this paper, we present the design of a novel copolymer consisting of EB and ES monomers with varying charge densities and different segmental arrangements. We present various intra- and intermolecular structural properties of the PANI chains using all-atom molecular dynamics (MD) simulations. Herein, we present a detailed conformational free energy analysis to understand the conformational transitions of the PANI chains. Our results show increased radius of gyration (Rg) values with increased charge density. Furthermore, we also present the H-bonding, free energy analysis, reduced density gradient (RDG), and solvent-accessible surface area (SASA) values for the observed conformational transitions of PANI. Therefore, these observations are crucial in understanding the complex behavior of chains for designing target-specific polymeric materials.

Treatment of flexibility of protein backbone in simulations of protein-ligand interactions using steered molecular dynamics

To ensure that an external force can break the interaction between a protein and a ligand, the steered molecular dynamics simulation requires a harmonic restrained potential applied to the protein backbone. A usual practice is that all or a certain number of protein's heavy atoms or Cα atoms are fixed, being restrained by a small force. This present study reveals that while fixing both either all heavy atoms and or all Cα atoms is not a good approach, while fixing a too small number of few atoms sometimes cannot prevent the protein from rotating under the influence of the bulk water layer, and the pulled molecule may smack into the wall of the active site. We found that restraining the Cα atoms under certain conditions is more relevant. Thus, we would propose an alternative solution in which only the Cα atoms of the protein at a distance larger than 1.2 nm from the ligand are restrained. A more flexible, but not too flexible, protein will be expected to lead to a more natural release of the ligand.

Improving the Potential of Mean Force and Nonequilibrium Pulling Simulations by Simultaneous Alchemical Modifications

We present an approach combining alchemical modifications and physical-pathway methods to calculate absolute binding free energies. The employed physical-pathway method is either a stratified umbrella sampling to calculate a potential of mean force or nonequilibrium pulling. We devised two basic approaches: the simultaneous approach (S-approach), where, along the physical unbinding pathway, an alchemical transformation of ligand-protein interactions is installed and deinstalled, and the prior-plus-simultaneous approach (PPS-approach), where, prior to the physical-pathway simulation, an alchemical transformation of ligand-protein interactions is installed in the binding site and deinstalled during the physical-pathway simulation. Using a mutant of T4 lysozyme with a benzene ligand as an example, we show that installation and deinstallation of soft-core interactions concurrent with physical ligand unbinding (S-approach) allow successful potential of mean force calculations and nonequilibrium pulling simulations despite the problems posed by the occluded nature of the lysozyme binding pocket. Good agreement between the potential of the mean-force-based S-approach and double decoupling simulations as well as a remarkable efficiency and accuracy of the nonequilibrium-pulling-based S-approach is found. The latter turned out to be more compute-efficient than the potential of mean force calculation by approximately 70%. Furthermore, we illustrate the merits of reducing ligand-protein interactions prior to potential of mean force calculations using the murine double minute homologue protein MDM2 with a p53-derived peptide ligand (PPS-approach). Here, the problem of breaking strong interactions in the binding pocket is transferred to a prior alchemical transformation that reduces the free-energy barrier between the bound and unbound state in the potential of mean force. Besides, disentangling physical ligand displacement from the deinstallation of ligand-protein interactions was seen to allow a more uniform sampling of distance histograms in the umbrella sampling. In the future, physical ligand unbinding combined with simultaneous alchemical modifications may prove useful in the calculation of protein-protein binding free energies, where sampling problems posed by multiple, possibly sticky interactions and potential steric clashes can thus be reduced.

The prediction of protein-ligand unbinding for modern drug discovery

INTRODUCTION

Drug-target thermodynamic and kinetic information have perennially important roles in drug design. The prediction of protein-ligand unbinding, which can provide important kinetic information, in experiments continues to face great challenges. Uncovering protein-ligand unbinding through molecular dynamics simulations has become efficient and inexpensive with the progress and enhancement of computing power and sampling methods.

AREAS COVERED

In this review, various sampling methods for protein-ligand unbinding and their basic principles are firstly briefly introduced. Then, their applications in predicting aspects of protein-ligand unbinding, including unbinding pathways, dissociation rate constants, residence time and binding affinity, are discussed.

EXPERT OPINION

Although various sampling methods have been successfully applied in numerous systems, they still have shortcomings and deficiencies. Most enhanced sampling methods require researchers to possess a wealth of prior knowledge of collective variables or reaction coordinates. In addition, most systems studied at present are relatively simple, and the study of complex systems in real drug research remains greatly challenging. Through the combination of machine learning and enhanced sampling methods, prediction accuracy can be further improved, and some problems encountered in complex systems also may be solved.

α-Helical Antimicrobial Peptide Encapsulation and Release from Boron Nitride Nanotubes: A Computational Study

Introduction

Antimicrobial peptides are potential therapeutics as anti-bacteria, anti-viruses, anti-fungi, or anticancers. However, they suffer from a short half-life and drug resistance which limit their long-term clinical usage.

Methods

Herein, we captured the encapsulation of antimicrobial peptide HA-FD-13 into boron nitride nanotube (BNNT) (20,20) and its release due to subsequent insertion of BNNT (14,14) with molecular dynamics simulation.

Results

The peptide-BNNT (20,20) van der Waals (vdW) interaction energy decreased to −270 kcal·mol⁻¹ at the end of the simulation (15 ns). However, during the period of 0.2–1.8 ns, when half of the peptide was inside the nanotube, the encapsulation was paused due to an energy barrier in the vicinity of BNNT and subsequently the external intervention, such that the self-adjustment of the peptide allowed full insertion. The free energy of the encapsulation process was −200.12 kcal·mol⁻¹, suggesting that the insertion procedure occurred spontaneously.

Discussion

Once the BNNT (14,14) entered into the BNNT (20,20), the peptide was completely released after 83.8 ps. This revealed that the vdW interaction between the BNNT (14,14) and BNNT (20,20) was stronger than between BNNT (20,20) and the peptide; therefore, the BNNT (14,14) could act as a piston pushing the peptide outside the BNNT (20,20). Moreover, the sudden drop in the vdW energy between nanotubes to the value of the −1300 Kcal·mol⁻¹ confirmed the self-insertion of the BNNT (14,14) into the BNNT (20,20) and correspondingly the release of the peptide.

Nicotinic receptor pharmacology in silico: Insights and challenges

Nicotinic acetylcholine receptors (nAChR) are homo- or hetero-pentameric ligand-gated ion channels of the Cys-loop superfamily and play important roles in the nervous system and muscles. Studies on nAChR benefit from in silico modeling due to the lack of high-resolution structures for most receptor subtypes and challenges in experiments addressing the complex mechanism of activation involving allosteric sites. Although there is myriad of computational modeling studies on nAChR, the multitude of the methods and parameters used in these studies makes modeling nAChR a daunting task, particularly for the non-experts in the field. To address this problem, the modeling literature on Torpedo nAChR and α7 nAChR were focused on as examples of heteromeric and homomeric nAChR, and the key in silico modeling studies between the years 1995-2019 were concisely reviewed. This was followed by a critical analysis of these studies by comparing the findings with each other and with the emerging experimental and computational data on nAChR. Based on these critical analyses, suggestions were made to guide the future researchers in the field of in silico modeling of nAChR. This article is part of the special issue on 'Contemporary Advances in Nicotine Neuropharmacology'.

Molecular dynamics simulations provide molecular insights into the role of HLA-B51 in Behçet's disease pathogenesis

Behçet's disease is an inflammatory disorder of unknown etiology. Genetic tendency has an important role in its pathogenesis, and HLA-B51, a class I MHC antigen, has been recognized as the strongest susceptibility factor for Behçet's disease. Despite the confirmation of the association of HLA-B51 with Behçet's disease in different populations, its pathogenic mechanisms remain elusive. HLA-B51 differs in only two amino acids from HLA-B52, other split antigen of HLA-B5, which is not associated with Behçet's disease. These two amino acids are located in the B pocket of the antigen-binding groove, which occupies the second amino acids of the bound peptides. To understand the nature of the HLA-peptide interactions, differences in structure and dynamics of two HLA alleles were investigated by molecular dynamics simulations using YAYDGKDYI, LPRSTVINI, and IPYQDLPHL peptides. For HLA-B51, all bound peptides fluctuated to larger extent than HLA-B52. Free energy profiles of unbinding process for YAYDGKDYI by steered molecular dynamics simulations showed that unbinding from HLA-B52 results in greater free energy differences than HLA-B51. These results suggest the possibility of an instability of HLA-B51 associated with the repertoire of peptides, and this finding may provide significant insight to its pathogenic role in Behçet's disease.

In Silico Modeling of the α7 Nicotinic Acetylcholine Receptor: New Pharmacological Challenges Associated with Multiple Modes of Signaling

The α7 nicotinic acetylcholine receptor is a homopentameric ion-channel of the Cys-loop superfamily characterized by its low probability of opening, high calcium permeability, and rapid desensitization. The α7 receptor has been targeted for the treatment of the cognitive symptoms of schizophrenia, depression, and Alzheimer's disease, but it is also involved in inflammatory modulation as a part of the cholinergic anti-inflammatory pathway. Despite its functional importance, in silico studies of the α7 receptor cannot produce a general model explaining the structural features of receptor activation, nor predict the mode of action for various ligand classes. Two particular problems in modeling the α7 nAChR are the absence of a high-resolution structure and the presence of five potentially nonequivalent orthosteric ligand binding sites. There is wide variability regarding the templates used for homology modeling, types of ligands investigated, simulation methods, and simulation times. However, a systematic survey focusing on the methodological similarities and differences in modeling α7 has not been done. In this work, we make a critical analysis of the modeling literature of α7 nAChR by comparing the findings of computational studies with each other and with experimental studies under the main topics of structural studies, ligand binding studies, and comparisons with other nAChR. In light of our findings, we also summarize current problems in the field and make suggestions for future studies concerning modeling of the α7 receptor.

Catalytic transport of molecular cargo using diffusive binding along a polymer track

Transport at the molecular scale is a prerequisite for the development of future molecular factories. Here, we have designed oligoanionic molecular sliders on polycationic tracks that exploit Brownian motion and diffusive binding to transport cargo without using a chemical fuel. The presence of the polymer tracks increases the rate of bimolecular reactions between modified sliders by over two orders of magnitude. Molecular dynamics simulations showed that the sliders not only diffuse, but also jump and hop surprisingly efficiently along polymer tracks. Inspired by acetyl-coenzyme A transporting and delivering acetyl groups in many essential biochemical processes, we developed a new and unconventional type of catalytic transport involving sliders (including coenzyme A) picking up, transporting and selectively delivering molecular cargo. Furthermore, we show that the concept of diffusive binding can also be utilized for the spatially controlled transport of chemical groups across gels. This work represents a new concept for designing functional nanosystems based on random Brownian motion.

Detailed potential of mean force studies on host-guest systems from the SAMPL6 challenge

Accurately predicting receptor-ligand binding free energies is one of the holy grails of computational chemistry with many applications in chemistry and biology. Many successes have been reported, but issues relating to sampling and force field accuracy remain significant issues affecting our ability to reliably calculate binding free energies. In order to explore these issues in more detail we have examined a series of small host-guest complexes from the SAMPL6 blind challenge, namely octa-acids (OAs)-guest complexes and Curcurbit[8]uril (CB8)-guest complexes. Specifically, potential of mean force studies using umbrella sampling combined with the weighted histogram method were carried out on both systems with both known and unknown binding affinities. We find that using standard force fields and straightforward simulation protocols we are able to obtain satisfactory results, but that simply scaling our results allows us to significantly improve our predictive ability for the unknown test sets: the overall RMSD of the binding free energy versus experiment is reduced from 5.59 to 2.36 kcal/mol; for the CB8 test system, the RMSD goes from 8.04 to 3.51 kcal/mol, while for the OAs test system, the RSMD goes from 2.89 to 0.95 kcal/mol. The scaling approach was inspired by studies on structurally related known benchmark sets: by simply scaling, the RMSD was reduced from 6.23 to 1.19 kcal/mol and from 2.96 to 0.62 kcal/mol for the CB8 benchmark system and the OA benchmark system, respectively. We find this scaling procedure to correct absolute binding affinities to be highly effective especially when working across a "congeneric" series with similar charge states. It is less successful when applied to mixed ligands with varied charges and chemical characteristics, but improvement is still realized in the present case. This approach suggests that there are large systematic errors in absolute binding free energy calculations that can be straightforwardly accounted for using a scaling procedure. Random errors are still an issue, but near chemical accuracy can be obtained using the present strategy in select cases.

Theoretical and Experimental Studies on Inclusion Complexes of Pinostrobin and β-Cyclodextrins

Pinostrobin (PNS) belongs to the flavanone subclass of flavonoids which shows several biological activities such as anti-inflammatory, anti-cancerogenic, anti-viral and anti-oxidative effects. Similar to other flavonoids, PNS has a quite low water solubility. The purpose of this work is to improve the solubility and the biological activities of PNS by forming inclusion complexes with β-cyclodextrin (βCD) and its derivatives, heptakis-(2,6-di-O-methyl)-β-cyclodextrin (2,6-DMβCD) and (2-hydroxypropyl)-β-cyclodextrin (HPβCD). The A_L-type diagram of the phase solubility studies of PNS exhibited the formed inclusion complexes with the 1:1 molar ratio. Inclusion complexes were prepared by the freeze-drying method and were characterized by differential scanning calorimetry (DSC). Two-dimensional nuclear magnetic resonance (2D-NMR) and steered molecular dynamics (SMD) simulation revealed two different binding modes of PNS, i.e., its phenyl- (P-PNS) and chromone- (C-PNS) rings preferably inserted into the cavity of βCD derivatives whilst only one orientation of PNS, where the C-PNS ring is inside the cavity, was detected in the case of the parental βCD. All PNS/βCDs complexes had a higher dissolution rate than free PNS. Both PNS and its complexes significantly exerted a lowering effect on the IL-6 secretion in LPS-stimulated macrophages and showed a moderate cytotoxic effect against MCF-7 and HeLa cancer cell lines in vitro.

Opposing Intermolecular Tuning of Ca2+ Affinity for Calmodulin by Neurogranin and CaMKII Peptides

We investigated the impact of bound calmodulin (CaM)-target compound structure on the affinity of calcium (Ca²⁺) by integrating coarse-grained models and all-atomistic simulations with nonequilibrium physics. We focused on binding between CaM and two specific targets, Ca²⁺/CaM-dependent protein kinase II (CaMKII) and neurogranin (Ng), as they both regulate CaM-dependent Ca²⁺ signaling pathways in neurons. It was shown experimentally that Ca²⁺/CaM (holoCaM) binds to the CaMKII peptide with overwhelmingly higher affinity than Ca²⁺-free CaM (apoCaM); the binding of CaMKII peptide to CaM in return increases the Ca²⁺ affinity for CaM. However, this reciprocal relation was not observed in the Ng peptide (Ng_13–49), which binds to apoCaM or holoCaM with binding affinities of the same order of magnitude. Unlike the holoCaM-CaMKII peptide, whose structure can be determined by crystallography, the structural description of the apoCaM-Ng_13–49 is unknown due to low binding affinity, therefore we computationally generated an ensemble of apoCaM-Ng_13–49 structures by matching the changes in the chemical shifts of CaM upon Ng_13–49 binding from nuclear magnetic resonance experiments. Next, we computed the changes in Ca²⁺ affinity for CaM with and without binding targets in atomistic models using Jarzynski’s equality. We discovered the molecular underpinnings of lowered affinity of Ca²⁺ for CaM in the presence of Ng_13–49 by showing that the N-terminal acidic region of Ng peptide pries open the β-sheet structure between the Ca²⁺ binding loops particularly at C-domain of CaM, enabling Ca²⁺ release. In contrast, CaMKII peptide increases Ca²⁺ affinity for the C-domain of CaM by stabilizing the two Ca²⁺ binding loops. We speculate that the distinctive structural difference in the bound complexes of apoCaM-Ng_13–49 and holoCaM-CaMKII delineates the importance of CaM’s progressive mechanism of target binding on its Ca²⁺ binding affinities.

Tinker-OpenMM: Absolute and relative alchemical free energies using AMOEBA on GPUs

The capabilities of the polarizable force fields for alchemical free energy calculations have been limited by the high computational cost and complexity of the underlying potential energy functions. In this work, we present a GPU-based general alchemical free energy simulation platform for polarizable potential AMOEBA. Tinker-OpenMM, the OpenMM implementation of the AMOEBA simulation engine has been modified to enable both absolute and relative alchemical simulations on GPUs, which leads to a ∼200-fold improvement in simulation speed over a single CPU core. We show that free energy values calculated using this platform agree with the results of Tinker simulations for the hydration of organic compounds and binding of host-guest systems within the statistical errors. In addition to absolute binding, we designed a relative alchemical approach for computing relative binding affinities of ligands to the same host, where a special path was applied to avoid numerical instability due to polarization between the different ligands that bind to the same site. This scheme is general and does not require ligands to have similar scaffolds. We show that relative hydration and binding free energy calculated using this approach match those computed from the absolute free energy approach. © 2017 Wiley Periodicals, Inc.

Ligand diffusion in proteins via enhanced sampling in molecular dynamics

Computational simulations in biophysics describe the dynamics and functions of biological macromolecules at the atomic level. Among motions particularly important for life are the transport processes in heterogeneous media. The process of ligand diffusion inside proteins is an example of a complex rare event that can be modeled using molecular dynamics simulations. The study of physical interactions between a ligand and its biological target is of paramount importance for the design of novel drugs and enzymes. Unfortunately, the process of ligand diffusion is difficult to study experimentally. The need for identifying the ligand egress pathways and understanding how ligands migrate through protein tunnels has spurred the development of several methodological approaches to this problem. The complex topology of protein channels and the transient nature of the ligand passage pose difficulties in the modeling of the ligand entry/escape pathways by canonical molecular dynamics simulations. In this review, we report a methodology involving a reconstruction of the ligand diffusion reaction coordinates and the free-energy profiles along these reaction coordinates using enhanced sampling of conformational space. We illustrate the above methods on several ligand-protein systems, including cytochromes and G-protein-coupled receptors. The methods are general and may be adopted to other transport processes in living matter.

Key Residues and Phosphate Release Routes in the Saccharomyces cerevisiae Pho84 Transceptor: THE ROLE OF TYR179 IN FUNCTIONAL REGULATION

Pho84, a major facilitator superfamily (MFS) protein, is the main high-affinity P_i transceptor in Saccharomyces cerevisiae Although transport mechanisms have been suggested for other MFS members, the key residues and molecular events driving transport by P_i:H⁺ symporters are unclear. The current Pho84 transport model is based on the inward-facing occluded crystal structure of the Pho84 homologue PiPT in the fungus Piriformospora indica However, this model is limited by the lack of experimental data on the regulatory residues for each stage of the transport cycle. In this study, an open, inward-facing conformation of Pho84 was used to study the release of P_i A comparison of this conformation with the model for P_i release in PiPT revealed that Tyr¹⁷⁹ in Pho84 (Tyr¹⁵⁰ in PiPT) is not part of the P_i binding site. This difference may be due to a lack of detailed information on the P_i release step in PiPT. Molecular dynamics simulations of Pho84 in which a residue adjacent to Tyr¹⁷⁹, Asp¹⁷⁸, is protonated revealed a conformational change in Pho84 from an open, inward-facing state to an occluded state. Tyr¹⁷⁹ then became part of the binding site as was observed in the PiPT crystal structure. The importance of Tyr¹⁷⁹ in regulating P_i release was supported by site-directed mutagenesis and transport assays. Using trehalase activity measurements, we demonstrated that the release of P_i is a critical step for transceptor signaling. Our results add to previous studies on PiPT, creating a more complete picture of the proton-coupled P_i transport cycle of a transceptor.

Evaluation of protein-ligand affinity prediction using steered molecular dynamics simulations

In computational drug design, ranking a series of compound analogs in a manner that is consistent with experimental affinities remains a challenge. In this study, we evaluated the prediction of protein-ligand binding affinities using steered molecular dynamics simulations. First, we investigated the appropriate conditions for accurate predictions in these simulations. A conic harmonic restraint was applied to the system for efficient sampling of work values on the ligand unbinding pathway. We found that pulling velocity significantly influenced affinity predictions, but that the number of collectable trajectories was less influential. We identified the appropriate pulling velocity and collectable trajectories for binding affinity predictions as 1.25 Å/ns and 100, respectively, and these parameters were used to evaluate three target proteins (FK506 binding protein, trypsin, and cyclin-dependent kinase 2). For these proteins using our parameters, the accuracy of affinity prediction was higher and more stable when Jarzynski's equality was employed compared with the second-order cumulant expansion equation of Jarzynski's equality. Our results showed that steered molecular dynamics simulations are effective for predicting the rank order of ligands; thus, they are a potential tool for compound selection in hit-to-lead and lead optimization processes.

A Self-Adaptive Steered Molecular Dynamics Method Based on Minimization of Stretching Force Reveals the Binding Affinity of Protein-Ligand Complexes

Binding affinity prediction of protein–ligand complexes has attracted widespread interest. In this study, a self-adaptive steered molecular dynamics (SMD) method is proposed to reveal the binding affinity of protein–ligand complexes. The SMD method is executed through adjusting pulling direction to find an optimum trajectory of ligand dissociation, which is realized by minimizing the stretching force automatically. The SMD method is then used to simulate the dissociations of 19 common protein–ligand complexes which are derived from two homology families, and the binding free energy values are gained through experimental techniques. Results show that the proposed SMD method follows a different dissociation pathway with lower a rupture force and energy barrier when compared with the conventional SMD method, and further analysis indicates the rupture forces of the complexes in the same protein family correlate well with their binding free energy, which reveals the possibility of using the proposed SMD method to identify the active ligand.

Calculating binding free energies for protein-carbohydrate complexes

A variety of computational techniques may be applied to compute theoretical binding free energies for protein-carbohydrate complexes. Elucidation of the intermolecular interactions, as well as the thermodynamic effects, that contribute to the relative strength of receptor binding can shed light on biomolecular recognition, and the resulting initiation or inhibition of a biological process. Three types of free energy methods are discussed here, including MM-PB/GBSA, thermodynamic integration, and a non-equilibrium alternative utilizing SMD. Throughout this chapter, the well-known concanavalin A lectin is employed as a model system to demonstrate the application of these methods to the special case of carbohydrate binding.

Mechanism of amyloid-β fibril elongation

Amyloid-β is an intrinsically disordered protein that forms fibrils in the brains of patients with Alzheimer's disease. To explore factors that affect the process of fibril growth, we computed the free energy associated with disordered amyloid-β monomers being added to growing amyloid fibrils using extensive molecular dynamics simulations coupled with umbrella sampling. We find that the mechanisms of Aβ40 and Aβ42 fibril elongation have many features in common, including the formation of an obligate on-pathway β-hairpin intermediate that hydrogen bonds to the fibril core. In addition, our data lead to new hypotheses for how fibrils may serve as secondary nucleation sites that can catalyze the formation of soluble oligomers, a finding in agreement with recent experimental observations. These data provide a detailed mechanistic description of amyloid-β fibril elongation and a structural link between the disordered free monomer and the growth of amyloid fibrils and soluble oligomers.

Conformational changes in acetylcholine binding protein investigated by temperature accelerated molecular dynamics

Despite the large number of studies available on nicotinic acetylcholine receptors, a complete account of the mechanistic aspects of their gating transition in response to ligand binding still remains elusive. As a first step toward dissecting the transition mechanism by accelerated sampling techniques, we study the ligand-induced conformational changes of the acetylcholine binding protein (AChBP), a widely accepted model for the full receptor extracellular domain. Using unbiased Molecular Dynamics (MD) and Temperature Accelerated Molecular Dynamics (TAMD) simulations we investigate the AChBP transition between the apo and the agonist-bound state. In long standard MD simulations, both conformations of the native protein are stable, while the agonist-bound structure evolves toward the apo one if the orientation of few key sidechains in the orthosteric cavity is modified. Conversely, TAMD simulations initiated from the native conformations are able to produce the spontaneous transition. With respect to the modified conformations, TAMD accelerates the transition by at least a factor 10. The analysis of some specific residue-residue interactions points out that the transition mechanism is based on the disruption/formation of few key hydrogen bonds. Finally, while early events of ligand dissociation are observed already in standard MD, TAMD accelerates the ligand detachment and, at the highest TAMD effective temperature, it is able to produce a complete dissociation path in one AChBP subunit.

Computational and experimental insights into the mechanism of substrate recognition and feedback inhibition of protoporphyrinogen oxidase

Protoporphyrinogen IX oxidase (PPO; EC 1.3.3.4) is an essential enzyme catalyzing the last common step in the pathway leading to heme and chlorophyll biosynthesis. Great interest in PPO inhibitors arises from both its significance to agriculture and medicine. However, the discovery of PPO inhibitors with ultrahigh potency and selectivity is hampered due to lack of structural and mechanistic understanding about the substrate recognition, which remains a longstanding question central in porphyrin biology. To understand the mechanism, a novel binding model of protogen (protoporphyrinogen IX, the substrate) was developed through extensive computational simulations. Subsequently, amino acid residues that are critical for protogen binding identified by computational simulations were substituted by mutagenesis. Kinetic analyses of these mutants indicated that these residues were critical for protogen binding. In addition, the calculated free energies of protogen binding with these mutants correlated well with the experimental data, indicating the reasonability of the binding model. On the basis of this novel model, the fundamental mechanism of substrate recognition was investigated by performing potential of mean force (PMF) calculations, which provided an atomic level description of conformational changes and pathway intermediates. The free energy profile revealed a feedback inhibition mechanism of proto (protoporphyrin IX, the product), which was also in agreement with experimental evidence. The novel mechanistic insights obtained from this study present a new starting point for future rational design of more efficient PPO inhibitors based on the product-bound PPO structure.

Computational gibberellin-binding channel discovery unraveling the unexpected perception mechanism of hormone signal by gibberellin receptor

Gibberellins (GAs) are phytohormones essential for many developmental processes in plants. In this work, fundamental mechanism of hormone perception by receptor GID1 has been studied by performing computational simulations, revealing a new GA-binding channel of GID1 and a novel hormone perception mechanism involving only one conformational state of GID1. The novel hormone perception mechanism demonstrated here is remarkably different from the previously proposed/speculated mechanism [Murase et al., Nature 2008, 456, 459] involving two conformational states ("OPEN" and "CLOSED") of GID1. According to the new perception mechanism, GA acts as a "conformational stabilizer," rather than the previously speculated "allosteric inducer," to induce the recognition of protein DELLA by GID1. The novel mechanistic insights obtained in this study provide a new starting point for further studies on the detailed molecular mechanisms of GID1 interacting with DELLA and various hormones and for mechanism-based rational design of novel, potent growth regulators that target crops and ornamental plants.

Computational design of new Peptide inhibitors for amyloid Beta (aβ) aggregation in Alzheimer's disease: application of a novel methodology

Alzheimer's disease is the most common form of dementia. It is a neurodegenerative and incurable disease that is associated with the tight packing of amyloid fibrils. This packing is facilitated by the compatibility of the ridges and grooves on the amyloid surface. The GxMxG motif is the major factor creating the compatibility between two amyloid surfaces, making it an important target for the design of amyloid aggregation inhibitors. In this study, a peptide, experimentally proven to bind Aβ40 fibrils at the GxMxG motif, was mutated by a novel methodology that systematically replaces amino acids with residues that share similar chemical characteristics and subsequently assesses the energetic favorability of these mutations by docking. Successive mutations are combined and reassessed via docking to a desired level of refinement. This methodology is both fast and efficient in providing potential inhibitors. Its efficiency lies in the fact that it does not perform all possible combinations of mutations, therefore decreasing the computational time drastically. The binding free energies of the experimentally studied reference peptide and its three top scoring derivatives were evaluated as a final assessment/valuation. The potential of mean forces (PMFs) were calculated by applying the Jarzynski‚s equality to results of steered molecular dynamics simulations. For all of the top scoring derivatives, the PMFs showed higher binding free energies than the reference peptide substantiating the usage of the introduced methodology to drug design.

Toward quantitative estimates of binding affinities for protein-ligand systems involving large inhibitor compounds: a steered molecular dynamics simulation route

Understanding binding mechanisms between enzymes and potential inhibitors and quantifying protein-ligand affinities in terms of binding free energy is of primary importance in drug design studies. In this respect, several approaches based on molecular dynamics simulations, often combined with docking techniques, have been exploited to investigate the physicochemical properties of complexes of pharmaceutical interest. Even if the geometric properties of a modeled protein-ligand complex can be well predicted by computational methods, it is still challenging to rank with chemical accuracy a series of ligand analogues in a consistent way. In this article, we face this issue calculating relative binding free energies of a focal adhesion kinase, an important target for the development of anticancer drugs, with pyrrolopyrimidine-based ligands having different inhibitory power. To this aim, we employ steered molecular dynamics simulations combined with nonequilibrium work theorems for free energy calculations. This technique proves very powerful when a series of ligand analogues is considered, allowing one to tackle estimation of protein-ligand relative binding free energies in a reasonable time. In our cases, the calculated binding affinities are comparable with those recovered from experiments by exploiting the Michaelis-Menten mechanism with a competitive inhibitor.

Ligand Binding Pathway Elucidation for Cryptophane Host-Guest Complexes

Modeling binding pathways can provide insight into molecular recognition, including kinetic mechanisms, barriers to binding, and gating effects. This work represents a novel computational approach, Hopping Minima, for the determination of conformational transitions of single molecules as well as binding pathways for molecular complexes. The method begins by thoroughly sampling a set of conformational minima for a molecular system. The natural motions of the system are modeled using the normal modes of the sampled minima. The natural motions are utilized to connect conformational minima and are finally combined to form association/binding pathways in the case of molecular complexes. We provide an implementation and example application of the method using alanine dipeptide and a set of chemical host-guest systems: two cryptophane hosts with two guest cations, trimethylammonium and tetramethylammonium. Our results demonstrate that conformational transitions can be modeled and extended to find binding pathways as well as energetic information relevant to the minimum conformations involved. This approach has advantages over simulation-based methods for studying systems with slow binding processes and can help design molecules with preferred binding kinetics.

Delineation of the unbinding pathway of α-conotoxin ImI from the α7 nicotinic acetylcholine receptor

α-Conotoxins potently and specifically inhibit isoforms of nicotinic acetylcholine receptors (nAChRs) and are used as molecular probes and as drugs or drug leads. Interactions occurring during binding and unbinding events are linked to binding kinetics, and knowledge of these interactions could help in the development of α-conotoxins as drugs. Here, the unbinding process for the prototypical α-conotoxin ImI/α7-nAChR system was investigated theoretically, and three exit routes were identified using random accelerated molecular dynamics simulations. The route involving the smallest conformation perturbation was further divided into three subpathways, which were studied using steered molecular dynamics simulations. Of the three subpathways, two had better experimental support and lower potential of mean force, indicating that they might be sampled more frequently. Additionally, these subpathways were supported by previous experimental studies. Several pairwise interactions, including a cation-π interaction and charge and hydrogen bond interactions, were identified as potentially playing important roles in the unbinding event.

Investigation of indazole unbinding pathways in CYP2E1 by molecular dynamics simulations

Human microsomal cytochrome P450 2E1 (CYP2E1) can oxidize not only low molecular weight xenobiotic compounds such as ethanol, but also many endogenous fatty acids. The crystal structure of CYP2E1 in complex with indazole reveals that the active site is deeply buried into the protein center. Thus, the unbinding pathways and associated unbinding mechanisms remain elusive. In this study, random acceleration molecular dynamics simulations combined with steered molecular dynamics and potential of mean force calculations were performed to identify the possible unbinding pathways in CYP2E1. The results show that channel 2c and 2a are most likely the unbinding channels of CYP2E1. The former channel is located between helices G and I and the B-C loop, and the latter resides between the region formed by the F-G loop, the B-C loop and the β1 sheet. Phe298 and Phe478 act as the gate keeper during indazole unbinding along channel 2c and 2a, respectively. Previous site-directed mutagenesis experiments also supported these findings.

Antigen-antibody interactions and structural flexibility of a femtomolar-affinity antibody

The femtomolar-affinity mutant antibody (4M5.3) generated by directed evolution is interesting because of the potential of antibody engineering. In this study, the mutant and its wild type (4-4-20) were compared in terms of antigen-antibody interactions and structural flexibility to elucidate the effects of directed evolution. For this purpose, multiple steered molecular dynamics (SMD) simulations were performed. The pulling forces of SMD simulations elucidated the regions that form strong attractive interactions in the binding pocket. Structural analysis in these regions showed two important mutations for improving attractive interactions. First, mutation of Tyr102(H) to Ser (sequence numbering of Protein Data Bank entry 1FLR ) played a role in resolving the steric hindrance on the pathway of the antigen in the binding pocket. Second, mutation of Asp31(H) to His played a role in resolving electrostatic repulsion. Potentials of mean force (PMFs) of both the wild type and the mutant showed landscapes that do not include obvious intermediate states and go directly to the bound state. These landscapes were regarded as funnel-like binding free energy landscapes. Furthermore, the structural flexibility based on the fluctuations of the positions of atoms was analyzed. It was shown that the fluctuations in the positions of the antigen and residues in contact with antigen tend to be smaller in the mutant than in the wild type. This result suggested that structural flexibility decreases as affinity is improved by directed evolution. This suggestion is similar to the relationship between affinity and flexibility for in vivo affinity maturation, which was suggested by Romesberg and co-workers [Jimenez, R., et al. (2003) Proc. Natl. Acad. Sci. U.S.A.100, 92-97]. Consequently, the relationship was found to be applicable up to femotomolar affinity levels.

Slow dissociation of a charged ligand: analysis of the primary quinone Q(A) site of photosynthetic bacterial reaction centers

Reaction centers (RCs) are integral membrane proteins that undergo a series of electron transfer reactions during the process of photosynthesis. In the Q(A) site of RCs from Rhodobacter sphaeroides, ubiquinone-10 is reduced, by a single electron transfer, to its semiquinone. The neutral quinone and anionic semiquinone have similar affinities, which is required for correct in situ reaction thermodynamics. A previous study showed that despite similar affinities, anionic quinones associate and dissociate from the Q(A) site at rates ≈10(4) times slower than neutral quinones indicating that anionic quinones encounter larger binding barriers (Madeo, J.; Gunner, M. R. Modeling binding kinetics at the Q(A) site in bacterial reaction centers. Biochemistry 2005, 44, 10994-11004). The present study investigates these barriers computationally, using steered molecular dynamics (SMD) to model the unbinding of neutral ground state ubiquinone (UQ) and its reduced anionic semiquinone (SQ(-)) from the Q(A) site. In agreement with experiment, the SMD unbinding barrier for SQ(-) is larger than for UQ. Multi Conformational Continuum Electrostatics (MCCE), used here to calculate the binding energy, shows that SQ(-) and UQ have comparable affinities. In the Q(A) site, there are stronger binding interactions for SQ(-) compared to UQ, especially electrostatic attraction to a bound non-heme Fe(2+). These interactions compensate for the higher SQ(-) desolvation penalty, allowing both redox states to have similar affinities. These additional interactions also increase the dissociation barrier for SQ(-) relative to UQ. Thus, the slower SQ(-) dissociation rate is a direct physical consequence of the additional binding interactions required to achieve a Q(A) site affinity similar to that of UQ. By a similar mechanism, the slower association rate is caused by stronger interactions between SQ(-) and the polar solvent. Thus, stronger interactions for both the unbound and bound states of charged and highly polar ligands can slow their binding kinetics without a conformational gate. Implications of this for other systems are discussed.

Charybdotoxin unbinding from the mKv1.3 potassium channel: a combined computational and experimental study

Charybdotoxin, belonging to the group of so-called scorpion toxins, is a short peptide able to block many voltage-gated potassium channels, such as mKv1.3, with high affinity. We use a reliable homology model based on the high-resolution crystal structure of the 94% sequence identical homologue Kv1.2 for charybdotoxin docking followed by molecular dynamics simulations to investigate the mechanism and energetics of unbinding, tracing the behavior of the channel protein and charybdotoxin during umbrella-sampling simulations as charybdotoxin is moved away from the binding site. The potential of mean force is constructed from the umbrella sampling simulations and combined with K(d) and free energy values gained experimentally using the patch-clamp technique to study the free energy of binding at different ion concentrations and the mechanism of the charybdotoxin-mKv1.3 binding process. A possible charybdotoxin binding mechanism is deduced that includes an initial hydrophobic contact followed by stepwise electrostatic interactions and finally optimization of hydrogen bonds and salt bridges.

Biophysical and computational studies of membrane penetration by the GRP1 pleckstrin homology domain

The pleckstrin homology (PH) domain of the general receptor for phosphoinositides 1 (GRP1) exhibits specific, high-affinity, reversible binding to phosphatidylinositol (3,4,5)-trisphosphate (PI(3,4,5)P(3)) at the plasma membrane, but the nature and extent of the interaction between this bound complex and the surrounding membrane environment remains unclear. Combining equilibrium and nonequilibrium molecular dynamics (MD) simulations, NMR spectroscopy, and monolayer penetration experiments, we characterize the membrane-associated state of GRP1-PH. MD simulations show loops flanking the binding site supplement the interaction with PI(3,4,5)P(3) through multiple contacts with the lipid bilayer. NMR data show large perturbations in chemical shift for these loop regions on binding to PI(3,4,5)P(3)-containing DPC micelles. Monolayer penetration experiments and further MD simulations demonstrate that mutating hydrophobic residues to polar residues in the flanking loops reduces membrane penetration. This supports a "dual-recognition" model of binding, with specific GRP1-PH-PI(3,4,5)P(3) interactions supplemented by interactions of loop regions with the lipid bilayer.

Accurate determination of the binding free energy for KcsA-charybdotoxin complex from the potential of mean force calculations with restraints

Free energy calculations for protein-ligand dissociation have been tested and validated for small ligands (50 atoms or less), but there has been a paucity of studies for larger, peptide-size ligands due to computational limitations. Previously we have studied the energetics of dissociation in a potassium channel-charybdotoxin complex by using umbrella sampling molecular-dynamics simulations, and established the need for carefully chosen coordinates and restraints to maintain the physiological ligand conformation. Here we address the ligand integrity problem further by constructing additional potential of mean forces for dissociation of charybdotoxin using restraints. We show that the large discrepancies in binding free energy arising from simulation artifacts can be avoided by using appropriate restraints on the ligand, which enables determination of the binding free energy within the chemical accuracy. We make several suggestions for optimal choices of harmonic potential parameters and restraints to be used in binding studies of large ligands.