Molecular dynamics simulations of bovine rhodopsin: influence of protonation states and different membrane-mimicking environments

G-protein coupled receptors (GPCRs) are a protein family of outstanding pharmaceutical interest. GPCR homology models, based on the crystal structure of bovine rhodopsin, have been shown to be valuable tools in the drug-design process. The initial model is often refined by molecular dynamics (MD) simulations, a procedure that has been recently discussed controversially. We therefore analyzed MD simulations of bovine rhodopsin in order to identify contacts that could serve as constraints in the simulation of homology models. Additionally, the effect of an N-terminal truncation, the nature of the membrane mimic, the influence of varying protonation states of buried residues and the importance of internal water molecules was analyzed. All simulations were carried out using the program-package GROMACS. While N-terminal truncation negatively influenced the overall protein stability, a stable simulation was possible in both solvent environments. As regards the protonation state of titratable sites, the experimental data could be reproduced by the program UHBD (University of Houston Brownian Dynamics), suggesting its application for studying homology models of GPCRs. A high flexibility was observed for internal water molecules at some sites. Finally, interhelical hydrogen-bonding interactions could be derived, which can now serve as constraints in the simulations of GPCR homology models.     

Calculation and Visualization of Atomistic Mechanical Stresses in Nanomaterials and Biomolecules

Many biomolecules have machine-like functions, and accordingly are discussed in terms of mechanical properties like force and motion. However, the concept of stress, a mechanical property that is of fundamental importance in the study of macroscopic mechanics, is not commonly applied in the biomolecular context. We anticipate that microscopical stress analyses of biomolecules and nanomaterials will provide useful mechanistic insights and help guide molecular design. To enable such applications, we have developed Calculator of Atomistic Mechanical Stress (CAMS), an open-source software package for computing atomic resolution stresses from molecular dynamics (MD) simulations. The software also enables decomposition of stress into contributions from bonded, nonbonded and Generalized Born potential terms. CAMS reads GROMACS topology and trajectory files, which are easily generated from AMBER files as well; and time-varying stresses may be animated and visualized in the VMD viewer. Here, we review relevant theory and present illustrative applications.     

Generalized ordering parameter S2 for N-H peptide bonding as measure of the conformational flexibility of proteins: comparison of algorithms of S2 calculation from molecular dynamics simulation data

NMR spectroscopy in combination with molecular dynamics (MD) simulation technique is a powerful tool for investigation of biomolecule conformations. In the present study we compare the two algorithms of calculation of the generalized order parameter S2 for the peptide N-H bond from MD simulation data. The program based on this algorithm has been written and the effect of an adjustible parameter on results of calculation was studied. HIV-1 protease MD simulation was performed during 7.36 ns using GROMACS 3.1.4. software package. The calculated order parameter values are in good agreement with those obtained by other authors using NMR.     

Gramicidin A channel as a test ground for molecular dynamics force fields

We use the well-known structural and functional properties of the gramicidin A channel to test the appropriateness of force fields commonly used in molecular dynamics (MD) simulations of ion channels. For this purpose, the high-resolution structure of the gramicidin A dimer is embedded in a dimyristoylphosphatidylcholine bilayer, and the potential of mean force of a K(+) ion is calculated along the channel axis using the umbrella sampling method. Calculations are performed using two of the most common force fields in MD simulations: CHARMM and GROMACS. Both force fields lead to large central barriers for K(+) ion permeation, that are substantially higher than those deduced from the physiological data by inverse methods. In long MD simulations lasting over 60 ns, several ions are observed to enter the binding site but none of them crossed the channel despite the presence of a large driving field. The present results, taken together with many earlier studies, highlights the shortcomings of the standard force fields used in MD simulations of ion channels and calls for construction of more appropriate force fields for this purpose.     

Molecular dynamics simulation of a decasaccharide fragment of heparin in aqueous solution

Molecular dynamics (MD) simulations on heparin-water-sodium systems were carried out in order to establish a simulation protocol able to represent heparin solution conformation under physiological conditions. Atomic charges suitable for heparin oligosaccharides were obtained from ab initio quantum-mechanical computations, at the 6-31G(**) level. The GROMACS forcefield, the SPC, and SPC/E water models were employed. Also heparin was simulated with IdoA residues in 1C(4) or 2S(0) conformational states. The results of the performed MD simulations are in agreement with the available experimental data, suggesting that this approach can be applied for the study of heparin interactions with its target proteins and thus play a role in the development of new antithrombotic agents.     

Experimental validation of molecular dynamics simulations of lipid bilayers: a new approach

A novel protocol has been developed for comparing the structural properties of lipid bilayers determined by simulation with those determined by diffraction experiments, which makes it possible to test critically the ability of molecular dynamics simulations to reproduce experimental data. This model-independent method consists of analyzing data from molecular dynamics bilayer simulations in the same way as experimental data by determining the structure factors of the system and, via Fourier reconstruction, the overall transbilayer scattering-density profiles. Multi-nanosecond molecular dynamics simulations of a dioleoylphosphatidylcholine bilayer at 66% RH (5.4 waters/lipid) were performed in the constant pressure and temperature ensemble using the united-atom GROMACS and the all-atom CHARMM22/27 force fields with the GROMACS and NAMD software packages, respectively. The quality of the simulated bilayer structures was evaluated by comparing simulation with experimental results for bilayer thickness, area/lipid, individual molecular-component distributions, continuous and discrete structure factors, and overall scattering-density profiles. Neither the GROMACS nor the CHARMM22/27 simulations reproduced experimental data within experimental error. The widths of the simulated terminal methyl distributions showed a particularly strong disagreement with the experimentally observed distributions. A comparison of the older CHARMM22 with the newer CHARMM27 force fields shows that significant progress is being made in the development of atomic force fields for describing lipid bilayer systems empirically.     

Virtual substitution scan via single‐step free energy perturbation

With the rapid expansion of our computing power, molecular dynamics (MD) simulations ranging from hundreds of nanoseconds to microseconds or even milliseconds have become increasingly common. The majority of these long trajectories are obtained from plain (vanilla) MD simulations, where no enhanced sampling or free energy calculation method is employed. To promote the “recycling” of these trajectories, we developed the Virtual Substitution Scan (VSS) toolkit as a plugin of the open‐source visualization and analysis software VMD. Based on the single‐step free energy perturbation (sFEP) method, VSS enables the user to post‐process a vanilla MD trajectory for a fast free energy scan of substituting aryl hydrogens by small functional groups. Dihedrals of the functional groups are sampled explicitly in VSS, which improves the performance of the calculation and is found particularly important for certain groups. As a proof‐of‐concept demonstration, we employ VSS to compute the solvation free energy change upon substituting the hydrogen of a benzene molecule by 12 small functional groups frequently considered in lead optimization. Additionally, VSS is used to compute the relative binding free energy of four selected ligands of the T4 lysozyme. Overall, the computational cost of VSS is only a fraction of the corresponding multi‐step FEP (mFEP) calculation, while its results agree reasonably well with those of mFEP, indicating that VSS offers a promising tool for rapid free energy scan of small functional group substitutions. © 2016 Wiley Periodicals, Inc. Biopolymers 105: 324–336, 2016.     

Docking and molecular dynamics simulations of peroxisome proliferator activated receptors interacting with pan agonist sodelglitazar

PPAR (peroxisome proliferator-activated receptor) pan agonists play a critical role in treating metabolic diseases, especially the Type-2 diabetes mellitus (T2DM). GlaxoSmithKline's sodelglitazar (GW677954) is one of the potent PPAR pan agonists, which is currently being investigated in Phase II clinical trials for the treatment of T2DM and its complications. The present study was aimed at investigation into the effect of sodelglitazar at the binding pockets of PPARs. The Schrodinger Suite program (2009) was used for the molecular docking, while the GROMACS program used for the molecular dynamics (MD) simulations. The results thus obtained showed that sodelglitazar being docked well in the active site of PPARs. It was revealed by the MD simulations that the structures of the receptors remained quite stable during the simulations and that the important AF-2 helix showed less flexibility after binding with sodelglitazar. Also, it was observed that sodelglitazar could periodically form hydrogen bonds with the AF-2 helix of PPARs to stabilize the AF-2 helix in an active conformation. Our findings have confirmed that GlaxoSmithKline's sodelglitazar can activate the PPARs, which is quite consistent with the previous biological studies.     

Study of the Classical Cassie Theory and Wenzel Theory Used in Nanoscale

Recently,there are hesitations in the application scope of the classical Cassie theory and Wenzel theory.In this paper,Molecular Dynamics(MD)simulations are use...     

FitDockApp：基于模板的分子对接PyMOL插件图形化计算软件

目的 分子对接在预测分子之间的结合模式和亲和力方面起着至关重要的作用,是计算结构生物学和计算机辅助药物设计研究的重要方法。本研究团队近期开发了一款基于模板的新型对接方法FitDock,当存在近似的蛋白质配体模板时,它在准确性和速度方面都超过了业界常用的分子对接方法。为了增强FitDock方法的可用性,使其在分子模拟领域得到更广泛的应用,很有必要发展图像化的软件工具。方法 基于Python图像化编程,本文开发了FitDockApp,这是分子可视化软件PyMOL的插件软件。结果 FitDockApp能够通过操作窗口界面,实现基于模板的分子对接和配体结构比对,实时显示预测三维结构,并提供将对接文件上传到实验室服务器获取最优模板的便利。此外,FitDockApp还具备批量对接功能。结论 FitDockApp通过用户友好的界面简化了对接过程,并提供丰富的功能,帮助研究人员获得精确的对接结果。FitDockApp是一款免费软件,兼容Windows和Linux系统,可在http://cao.labshare.cn/fitdock/下载。     

Theoretical studies on the interaction of biphenyl inhibitors with Mycobacterium tuberculosis protein tyrosine phosphatase MptpB

MptpB is an essential secreted virulence factor for M. tuberculosis. Inhibition of MptpB impairs mycobacterial survival in host macrophages and thus helps reduce tuberculosis infections. However, the binding mode of the biphenyl inhibitors, which are known as some of the most potent MptpB inhibitors, remains unclear. In this study, to understand the interactions between biphenyl inhibitors and MptpB, docking and molecular dynamics simulations were carried out using AutoDock and GROMACS softwares. Calculation results show that all the biphenyl inhibitors can be docked to the binding site of MptpB, with the acid warheads forming a hydrogen bond network at the active site. But the binding modes of other terminals of these inhibitors are different. The cyclohexyl and trifluoromethyl substituents at R1 and R2 sites are necessary for the inhibitors to adopt their double-site binding mechanism. The estimated binding affinities are basically consistent with the experimental results. MD simulations show that these binding complexes display different stability.     

Conformational transitions of flanking purines in HIV-1 RNA dimerization initiation site kissing complexes studied by CHARMM explicit solvent molecular dynamics

Dimerization of HIV-1 genomic RNA is initiated by kissing loop interactions at the Dimerization Initiation Site (DIS). Dynamics of purines that flank the 5' ends of the loop-loop helix in HIV-1 DIS kissing complex were explored using explicit solvent molecular dynamics (MD) simulations with the CHARMM force field. Multiple MD simulations (200 ns in total) of X-ray structures for HIV-1 DIS Subtypes A, B, and F revealed conformational variability of flanking purines. In particular, the flanking purines, which in the starting X-ray structures are bulged-out and stack in pairs, formed a consecutive stack of four bulged-out adenines at the beginning of several simulations. This conformation is seen in the crystal structure of DIS Subtype F with no interference from crystal packing, and was frequently reported in our preceding MD studies performed with the AMBER force field. However, as CHARMM simulations progressed, the four continuously stacked adenines showed conformational transitions from the bulged-out into the bulged-in geometries. Although such an arrangement has not been seen in any X-ray structure, it has been suggested by a recent NMR investigation. In CHARMM simulations, in the longer time scale, the flanking purines display the tendency to move to bulged-in conformations. This is in contrast with the AMBER simulations, which indicate a modest prevalence for bulged-out flanking base positions in line with the X-ray data. The simulations also suggest that the intermolecular stacking between purines from the opposite hairpins can additionally stabilize the kissing complex.     

Molecular dynamics investigations of the polysaccharide scleroglucan: first study on the triple helix structure

Explicit solvent molecular dynamics (MD) simulations on the triple helix of the polysaccharide Scleroglucan (Sclg) at two temperatures (273 and 300 K) were carried out. Owing to the complexity of the system, a united-atom force field, based on the properly modified GROMACS parameters, was adopted. To test these parameters for our system, MD simulations of the two disaccharidic units, representing the main chain and the side-chain linkages of the Sclg repeating unit, were performed and the results were compared with the literature data. The simulated triple helix of Sclg retained the main experimentally determined features of the polymer. The residence times of the solvent molecules at 273 and 300 K were analyzed. The results show that the more internal water molecules, interacting with the core of the Sclg triplex are not influenced substantially by changing the temperature, on the contrary the water molecules, interacting with the side-chain glucose residues show more significant differences. These data suggest that the more external water molecules, interacting with the side chain, play a major role in the conformational transition experimentally observed at low temperature.     

Conformational space of alpha 1-4 glycosidic linkage: a molecular dynamics study

Molecular Dynamics simulations have been carried out for 100 ps on crystal structure of beta-cyclodextrin in vacuo and with explicit inclusion of solvent at constant pressure and constant temperature using the GROMOS MD algorithm, with a time step of 0.005 ps. The conformational space of the glycosidic linkage was studied by calculating two virtual dihedrals connecting the successive glucose units for the 2000 structures saved during the two simulations. Three preferred regions for alpha 1-4 glycosidic linkage were found in both the simulations. The use of these virtual dihedral angles in representing the glycosidic linkage is also brought out from these studies.     

Combining optimal control theory and molecular dynamics for protein folding

A new method to develop low-energy folding routes for proteins is presented. The novel aspect of the proposed approach is the synergistic use of optimal control theory with Molecular Dynamics (MD). In the first step of the method, optimal control theory is employed to compute the force field and the optimal folding trajectory for the Cα atoms of a Coarse-Grained (CG) protein model. The solution of this CG optimization provides an harmonic approximation of the true potential energy surface around the native state. In the next step CG optimization guides the MD simulation by specifying the optimal target positions for the Cα atoms. In turn, MD simulation provides an all-atom conformation whose Cα positions match closely the reference target positions determined by CG optimization. This is accomplished by Targeted Molecular Dynamics (TMD) which uses a bias potential or harmonic restraint in addition to the usual MD potential. Folding is a dynamical process and as such residues make different contacts during the course of folding. Therefore CG optimization has to be reinitialized and repeated over time to accomodate these important changes. At each sampled folding time, the active contacts among the residues are recalculated based on the all-atom conformation obtained from MD. Using the new set of contacts, the CG potential is updated and the CG optimal trajectory for the Cα atoms is recomputed. This is followed by MD. Implementation of this repetitive CG optimization-MD simulation cycle generates the folding trajectory. Simulations on a model protein Villin demonstrate the utility of the method. Since the method is founded on the general tools of optimal control theory and MD without any restrictions, it is widely applicable to other systems. It can be easily implemented with available MD software packages.     

Hybrid approaches to molecular simulation

Molecular dynamics (MD) simulation is an established method for studying the conformational changes that are important for protein function. Recent advances in hardware and software have allowed MD simulations over the same timescales as experiment, improving the agreement between theory and experiment to a large extent. However, running such simulations are costly, in terms of resources, storage, and trajectory analysis. There is still a place for techniques that involve short MD simulations. In order to overcome the sampling paucity of short time-scales, hybrid methods that include some form of MD simulation can exploit certain features of the system of interest, often combining experimental information in surprising ways. Here, we review some recent hybrid approaches to the simulation of proteins.