Lightweight object oriented structure analysis: Tools for building tools to analyze molecular dynamics simulations

LOOS (Lightweight Object Oriented Structure‐analysis) is a C++ library designed to facilitate making novel tools for analyzing molecular dynamics simulations by abstracting out the repetitive tasks, allowing developers to focus on the scientifically relevant part of the problem. LOOS supports input using the native file formats of most common biomolecular simulation packages, including CHARMM, NAMD, Amber, Tinker, and Gromacs. A dynamic atom selection language based on the C expression syntax is included and is easily accessible to the tool‐writer. In addition, LOOS is bundled with over 140 prebuilt tools, including suites of tools for analyzing simulation convergence, three‐dimensional histograms, and elastic network models. Through modern C++ design, LOOS is both simple to develop with (requiring knowledge of only four core classes and a few utility functions) and is easily extensible. A python interface to the core classes is also provided, further facilitating tool development. © 2014 Wiley Periodicals, Inc.     

Pyrite: A blender plugin for visualizing molecular dynamics simulations using industry‐standard rendering techniques

Molecular dynamics (MD) simulations provide critical insights into many biological mechanisms. Programs such as VMD, Chimera, and PyMOL can produce impressive simulation visualizations, but they lack many advanced rendering algorithms common in the film and video‐game industries. In contrast, the modeling program Blender includes such algorithms but cannot import MD‐simulation data. MD trajectories often require many gigabytes of memory/disk space, complicating Blender import. We present Pyrite, a Blender plugin that overcomes these limitations. Pyrite allows researchers to visualize MD simulations within Blender, with full access to Blender's cutting‐edge rendering techniques. We expect Pyrite‐generated images to appeal to students and non‐specialists alike. A copy of the plugin is available at http://durrantlab.com/pyrite/ , released under the terms of the GNU General Public License Version 3. © 2017 Wiley Periodicals, Inc.     

Investigation of chiral recognition by molecular micelles with molecular dynamics simulations

Molecular dynamics simulations were used to characterize the binding of the chiral drugs chlorthalidone and lorazepam to the molecular micelle poly-(sodium undecyl-(L)-leucine-valine). The project’s goal was to characterize the nature of chiral recognition in capillary electrophoresis separations that use molecular micelles as the chiral selector. The shapes and charge distributions of the chiral molecules investigated, their orientations within the molecular micelle chiral binding pockets, and the formation of stereoselective intermolecular hydrogen bonds with the molecular micelle were all found to play key roles in determining where and how lorazepam and chlorthalidone enantiomers interacted with the molecular micelle.     

Combination of COSMOmic and molecular dynamics simulations for the calculation of membrane–water partition coefficients

The importance of membrane‐water partition coefficients led to the recent extension of the conductor‐like screening model for realistic solvation (COSMO‐RS) to micelles and biomembranes termed COSMOmic. Compared to COSMO‐RS, this new approach needs structural information to account for the anisotropy of colloidal systems. This information can be obtained from molecular dynamics (MD) simulations. In this work, we show that this combination of molecular methods can efficiently be used to predict partition coefficients with good agreement to experimental data and enables screening studies. However, there is a discrepancy between the amount of data generated by MD simulations and the structural information needed for COSMOmic. Therefore, a new scheme is presented to extract data from MD trajectories for COSMOmic calculations. In particular, we show how to calculate the system structure from MD, the influence of lipid conformers, the relation to the COSMOmic layer size, and the water/lipid ratio impact. For a 1,2‐dimyristoyl‐sn‐glycero‐3‐phosphocholine (DMPC) bilayer, 66 partition coefficients for various solutes were calculated. Further, 52 partition coefficients for a 1‐palmitoyl‐2‐oleoyl‐sn‐glycero‐3‐phosphocholine (POPC) bilayer system were calculated. All these calculations were compared to experimental data. © 2013 Wiley Periodicals, Inc.     

Force-field parametrization and molecular dynamics simulations of Congo red

Congo red, a diazo dye widely used in medical diagnosis, is known to form supramolecular systems in solution. Such a supramolecular system may interact with various proteins. In order to examine the nature of such complexes empirical force field parameters for the Congo red molecule were developed. The parametrization of bonding terms closely followed the methodology used in the development of the charmm22 force field, except for the calculation of charges. Point charges were calculated from a fit to a quantum mechanically derived electrostatic potential using the CHELP-BOW method. Obtained parameters were tested in a series of molecular dynamics simulations of both a single molecule and a micelle composed of Congo red molecules. It is shown that newly developed parameters define a stable minimum on the hypersurface of the potential energy and crystal and ab initio geometries and rotational barriers are well reproduced. Furthermore, rotations around C-N bonds are similar to torsional vibrations observed in crystals of diphenyl-diazene, which confirms that the flexibility of the molecule is correct. Comparison of results obtained from micelles molecular dynamics simulations with experimental data shows that the thermal dependence of micelle creation is well reproduced.     

A box-counting-based algorithm for computing Shannon entropy in molecular dynamics simulations

A box-counting-based algorithm (SEBC) has been developed for the numerical computation of the Shannon entropy from samples of continuous functions. Its performance was tested by applying it to several samples of known continuous distribution functions. The results obtained with SEBC reproduced those obtained by analytical or numerical integration. SEBC was also employed for computing the Shannon entropies of the steric energy, Sh(E(S)), of several amino acids from their in vacuo NVE molecular dynamics simulations using the AMBER-4 force field. The results obtained correlate linearly with the experimental standard thermodynamic entropies of these compounds. This work points to the possibility of introducing straightforward and reliable calculations of thermodynamic entropies from empirical linear relationships with Sh(E(S)) obtained from MD simulations.     

Probing molecular mechanism of inhibitor bindings to bromodomain-containing protein 4 based on molecular dynamics simulations and principal component analysis

ABSTRACT

It is well known that bromodomain-containing protein 4 (BRD4) has been thought as a promising target utilized for treating various human diseases, such as inflammatory disorders, malignant tumours, acute myelogenous leukaemia (AML), bone diseases, etc. For this study, molecular dynamics (MD) simulations, binding free energy calculations, and principal component analysis (PCA) were integrated together to uncover binding modes of inhibitors 8P9, 8PU, and 8PX to BRD4(1). The results obtained from binding free energy calculations show that van der Waals interactions act as the main regulator in bindings of inhibitors to BRD4(1). The information stemming from PCA reveals that inhibitor associations extremely affect conformational changes, internal dynamics, and movement patterns of BRD4(1). Residue-based free energy decomposition method was wielded to unveil contributions of independent residues to inhibitor bindings and the data signify that hydrogen bonding interactions and hydrophobic interactions are decisive factors affecting bindings of inhibitors to BRD4(1). Meanwhile, eight residues Trp81, Pro82, Val87, Leu92, Leu94, Cys136, Asn140, and Ile146 are recognized as the common hot interaction spots of three inhibitors with BRD4(1). The results from this work are expected to provide a meaningfully theoretical guidance for design and development of effective inhibitors inhibiting of the activity of BRD4.     

Bond detectors for molecular dynamics simulations,Part I: Hydrogen bonds

Charge sensitivity analysis in AMBER force‐field resolution has been used in quest for detectors of hydrogen bonds (HBs). The process of HB formation was investigated on ab initio classical trajectories (B3LYP/6‐31G*) of different nucleobase pairs. Several charge sensitivities, namely: electronegativity, hardness, Fukui function (FF), and polarization matrix, were analyzed. The global and constrained equilibria were considered. It was demonstrated that FF indices and polarization matrix elements are good detectors of HB formation. © 2013 Wiley Periodicals, Inc.     

Charge distributions for molecular dynamics simulations from self-consistent polarization method

Partial atomic charges are important force field parameters. They are usually computed by applying quantum-chemical calculations and the assumed population scheme. In this study polarization consistent scheme of deriving a charge distribution inside solute molecule is proposed. The environment effect is explicitly taken into account by distributing solvent molecules around the solute target. The performed analysis includes a few computational schemes (HF, MP2, B3LYP, and M026X), basis sets (cc-pvnz, n = 2, 3, …, 6), and electrostatically derived charge distributions (KS, CHELP, CHELPG, and HLY). It is demonstrated that the environment effect is very important and cannot be disregarded. The second solvation shell should be included to achieve the charge convergence. Huge corrections to charge distribution are due to induction and dispersion. The B3LYP/cc-pvqz level of theory is recommended for deriving the charges within self-consistent polarization scheme.     

Understanding the adsorption of branched polyamine on surface of gold nanoparticles by molecular dynamics simulations

Molecular dynamics simulations were conducted to characterize the adsorption behavior of branched polyethylenimine (br‐PEI) on the surface of gold nanoparticles (AuNPs). We observed the preferential adsorption of br‐PEI on the [111] surface of AuNPs. Furthermore, br‐PEI maintained a flat arrangement on the surface and wrapped the AuNPs in a stable manner, thereby blocking the adsorption of H₂O molecules and other free br‐PEI molecules. The model and computational results provide theoretical support for relevant experiments. Copyright © 2016 John Wiley & Sons, Ltd.     

Nonspecific interaction forces at water-membrane interface by forced molecular dynamics simulations

Nonspecific interactions are the main driving forces for the behavior of molecules with great affinity for biologic membranes. To investigate not only the molecular details of these interactions but to estimate their magnitude as well, the theoretical method of Forced Molecular Dynamics Simulations, based on the Atomic Force Spectroscopy experimental technique, was applied. In this approach, an additional one-dimensional elastic force, representing the cantilever probe, was incorporated to the force field of a Molecular Dynamics computational program. This force represents a spring fixed on one end to a selected atom of the molecule; the other end of the spring is displaced at constant velocity to pull the molecule out of the membrane. The force experimented by the molecule due to the spring, is proportional to the spring elongation relative to its equilibrium position. This value is registered during the entire simulation, and its maximum value will determine the molecule-membrane interaction force. Nonexplicit medium simulations were carried out. Polar and apolar media were considered according to their polarizability degree and a specific dielectric constant value was assigned. In this approach, the membrane was considered as the apolar region limited by two flat surfaces with a polar aqueous medium. The potential energy discontinuity at the interfaces was smoothed by considering the polarization-induced effects using the image method. The results of this methodology are presented using a small system, a single Alanine amino acid model, which enables extended simulations in a microsecond time scale. The confinement of this amino acid at the interface reduces its degrees of freedom and forces it to adopt one of the six defined conformations. A correlation between these stable structures at the water-membrane interface and the interaction force value was determined.     

Essential dynamics/factor analysis for the interpretation of molecular dynamics trajectories

Subject of this work is the analysis of molecular dynamics (MD) trajectories of neurophysins I (NPI) and II (NPII) and their complexes with the neurophyseal nonapeptide hormones oxytocin (OT) and vasopresssin (VP), respectively, simulated in water. NPs serve in the neurosecretory granules as carrier proteins for the hormones before their release to the blood. The starting data consisted of two pairs of different trajectories for each of the (NPII/VP)2 and (NPI/OT)2 heterotetramers and two more trajectories for the NPII2 and NPI2 homodimers (six trajectories in total). Using essential dynamics which, to our judgement, is equivalent to factor analysis, we found that only about 10 degrees of freedom per trajectory are necessary and sufficient to describe in full the motions relevant for the function of the protein. This is consistent with these motions to explain about 90% of the total variance of the system. These principal degrees of freedom represent slow anharmonic motional modes, clearly pointing at distinguished mobility of the atoms involved in the protein's functionality.     

A new set of molecular mechanics parameters for hydroxyproline and its use in molecular dynamics simulations of collagen-like peptides

Recently, the importance of proline ring pucker conformations in collagen has been suggested in the context of hydroxylation of prolines. The previous molecular mechanics parameters for hydroxyproline, however, do not reproduce the correct pucker preference. We have developed a new set of parameters that reproduces the correct pucker preference. Our molecular dynamics simulations of proline and hydroxyproline monomers as well as collagen-like peptides, using the new parameters, support the theory that the role of hydroxylation in collagen is to stabilize the triple helix by adjusting to the right pucker conformation (and thus the right phi angle) in the Y position.     

The distributed diagonal force decomposition method for parallelizing molecular dynamics simulations

Parallelization is an effective way to reduce the computational time needed for molecular dynamics simulations. We describe a new parallelization method, the distributed-diagonal force decomposition method, with which we extend and improve the existing force decomposition methods. Our new method requires less data communication during molecular dynamics simulations than replicated data and current force decomposition methods, increasing the parallel efficiency. It also dynamically load-balances the processors' computational load throughout the simulation. The method is readily implemented in existing molecular dynamics codes and it has been incorporated into the CHARMM program, allowing its immediate use in conjunction with the many molecular dynamics simulation techniques that are already present in the program. We also present the design of the Force Decomposition Machine, a cluster of personal computers and networks that is tailored to running molecular dynamics simulations using the distributed diagonal force decomposition method. The design is expandable and provides various degrees of fault resilience. This approach is easily adaptable to computers with Graphics Processing Units because it is independent of the processor type being used.     

A comparison of different initialization protocols to obtain statistically independent molecular dynamics simulations

We study how the results of molecular dynamics (MD) simulations are affected by various choices during the setup, e.g., the starting velocities, the solvation, the location of protons, the conformation of His, Asn, and Gln residues, the protonation and titration of His residues, and the treatment of alternative conformations. We estimate the binding affinity of ligands to four proteins calculated with the MM/GBSA method (molecular mechanics combined with a generalized Born and surface area solvation energy). For avidin and T4 lysozyme, all variations gave similar results within 2 kJ/mol. For factor Xa, differences in the solvation or in the selection of alternative conformations gave results that are significantly different from those of the other approaches by 4-6 kJ/mol, whereas for galectin-3, changes in the conformations, rotations, and protonation gave results that differed by 10 kJ/mol, but only if residues close to the binding site were modified. This shows that the results of MM/GBSA calculations are reasonably reproducible even if the MD simulations are set up with different software. Moreover, we show that the sampling of phase space can be enhanced by solvating the systems with different equilibrated water boxes, in addition to the common use of different starting velocities. If different conformations are available in the crystal structure, they should also be employed to enhance the sampling. Protonation, ionization, and conformations of Asn, Gln, and His may also be used to enhance sampling, but great effort should be spent to obtain as reliable predictions as possible close to the active site.