CHARMM: The biomolecular simulation program

CHARMM (Chemistry at HARvard Molecular Mechanics) is a highly versatile and widely used molecular simulation program. It has been developed over the last three decades with a primary focus on molecules of biological interest, including proteins, peptides, lipids, nucleic acids, carbohydrates, and small molecule ligands, as they occur in solution, crystals, and membrane environments. For the study of such systems, the program provides a large suite of computational tools that include numerous conformational and path sampling methods, free energy estimators, molecular minimization, dynamics, and analysis techniques, and model-building capabilities. The CHARMM program is applicable to problems involving a much broader class of many-particle systems. Calculations with CHARMM can be performed using a number of different energy functions and models, from mixed quantum mechanical-molecular mechanical force fields, to all-atom classical potential energy functions with explicit solvent and various boundary conditions, to implicit solvent and membrane models. The program has been ported to numerous platforms in both serial and parallel architectures. This article provides an overview of the program as it exists today with an emphasis on developments since the publication of the original CHARMM article in 1983. © 2009 Wiley Periodicals, Inc.J Comput Chem, 2009.Peer Reviewedhttp://deepblue.lib.umich.edu/bitstream/2027.42/63074/1/21287_ftp.pd

A mesoscale model for coarse-grained protein dynamics

Proteins are the essential units of biological processes, but modelling their dynamics is a very computationally expensive task. A wide variety of simulation techniques exist, a popular example being Molecular Dynamics. However, such models typically involve detailed simulation of the protein's structure at or near the atomic level and as such are unsuitable for modelling biological systems composed of large or multiply interacting proteins. This research takes a coarse-graining approach, called Fluctuating Finite Element Analysis, in which large, globular proteins are approximated by viscoelastic continua subject to thermal noise. Each protein is discretised into a tetrahedral mesh, parameterised locally by its bulk continuum properties. The forces are then calculated using Finite Element Analysis. A parallel implementation of the FFEA algorithm has been developed for use on high performance computing facilities. The scalability of the algorithm with respect to number of cores and system size, and its stability with respect to integration time step has been investigated. A pipeline for fully automated FFEA system creation from atomistic (X-ray crystallography and NMR) or low resolution data (cryo-EM and SAXS) has also been developed. In order to tackle multiprotein systems, the FFEA model has been extended to include van derWaals interactions and electrostatics. FFEA has been applied to a number of diverse biological systems. The van der Waals scheme was tested through simulation of myoglobins interacting with a polystyrene substrate. The major modes of motion of V- and A- type rotary ATPases were extracted using Principal Component Analysis, and compared with the normal modes obtained from the Elastic Network Model. Finally, the effect of axonemal dynein c's interaction with the microtubule track on its step length and exploration of binding sites was investigated. A mapping was developed to allow in-simulation conformational switching of the dynein motor

Multiscale simulations of protein dynamics

Characteristic timescales associated with the function of biomolecules, like proteins, range from femtoseconds up to minutes, whereas their corresponding spatial extent ranges from few ̊A to μm when associating in large macromolecular complexes. Moreover, biomolecules are functional in a large variety of different physico-chemical conditions strongly dependent on pH, ionic strength, crowding agents, etc. This huge complexity is hard to be studied with an arbitrary level of resolution embracing all these spatial and temporal scales. Molecular simulation is a well established approach to gain mechanistic insights into the function of biomolecular systems, producing atomistically detailed models of in vitro and/or in vivo conditions. I present in this thesis two projects that aim at improving on the current limitations of multiscale molecular simulations, namely (i) the sampling of large systems, and (ii) the detailed representation and description of realistic physiological conditions. Addressing the first issue, I propose a new coarse-grained model for proteins to be used in molecular dynamics simulations. This coarse-grained model is based on a more accurate description of protein electrostatics, which accounts for dipolar contributions. The parameterization of this force field is based on force-matching methods and on the use of a particle swarm optimization heuristic algorithm. The obtained results are encouraging being structural and electrostatic properties accurately reproduced with the coarse-grained model for a variety of protein folds. Moreover, the parameterization procedure can be straightforwardly applied to any protein, and can be extended to a larger dataset to generate a fully transferable coarse-grained force field, to be applied to any protein and any large macromolecular assemblies for which long all-atom simulations are still a challenge. While the development and use of coarse-grained models are important to tackle the limited sampling of large systems, it is still important to use all-atom molecular dynamics simulation to investigate with high accuracy in physiological conditions protein dynamics. For this reason, I present here the results of state-of-art molecular dynamics simulations applied to study the influence of crowding agents on the internal dynamics of the protein ubiquitin. Their analysis allows to describe how ubiquitin dynamics is slaved by crowding agents. This work demonstrates that the description of protein dynamics should take into account its intrinsic multiscale nature. The development and applications of coarse-grained models permit to simulate proteins at low computational cost, the use of atomistic simulations allows to accurately describe proteins in absence and presence of crowding agents, and both of them permit to highlight the essence of protein dynamics

Domain Decomposition Based Hybrid Methods of Finite Element and Finite Difference and Applications in Biomolecule Simulations

The dielectric continuum models, such as Poisson Boltzmann equation (PBE), size modified PBE (SMPBE), and nonlocal modified PBE (NMPBE), are important models in predicting the electrostatics of a biomolecule in an ionic solvent. To solve these dielectric continuum models efficiently, in this dissertation, new finite element and finite difference hybrid methods are constructed by Schwartz domain decomposition techniques based on a special seven-box partition of a cubic domain. As one important part of these methods, a finite difference optimal solver --- the preconditioned conjugate gradient method using a multigrid V-cycle preconditioner --- is described in details and proved to have a convergence rate independent of mesh size in solving a symmetric positive definite linear system. These new hybrid algorithms are programmed in Fortran, C, and Python based on the efficient finite element library DOLFIN from the FEniCS project, and are well validated by test models with known analytical solutions. Comparison numerical tests between the new hybrid solvers and the corresponding finite element solvers are done to show the improvement in efficiency. Finally, as applications, solvation free energy and binding free energy calculations are done and then compared to the experiment data

A Computational Model of Protein Induced Membrane Morphology with Geodesic Curvature Driven Protein-Membrane Interface

Continuum or hybrid modeling of bilayer membrane morphological dynamics induced by embedded proteins necessitates the identification of protein-membrane interfaces and coupling of deformations of two surfaces. In this article we developed (i) a minimal total geodesic curvature model to describe these interfaces, and (ii) a numerical one-one mapping between two surface through a conformal mapping of each surface to the common middle annulus. Our work provides the first computational tractable approach for determining the interfaces between bilayer and embedded proteins. The one-one mapping allows a convenient coupling of the morphology of two surfaces. We integrated these two new developments into the energetic model of protein-membrane interactions, and developed the full set of numerical methods for the coupled system. Numerical examples are presented to demonstrate (1) the efficiency and robustness of our methods in locating the curves with minimal total geodesic curvature on highly complicated protein surfaces, (2) the usefulness of these interfaces as interior boundaries for membrane deformation, and (3) the rich morphology of bilayer surfaces for different protein-membrane interfaces

Molecular modeling and simulation of pH effects in lipid bilayers

Tese de doutoramento, Bioquímica (Bioquímica Teórica), Universidade de Lisboa, Faculdade de Ciências, 2015Biological membranes are the first barrier between cell and the extracellular environment, and also their first via of interaction. These biological structures are complex and diverse but share a common feature: they are all supported by a lipid bilayer. This bilayer is often negative and sensitive to pH and ionic strength. Computational / theoretical methods have been used to understand how these two properties inuence the bilayer structure. However, the methods available at the beginning of our work had two significant bottlenecks that we attempted to get rid of. First, none of the available Poisson{Boltzmann (PB) solvers were able to deal with both periodicity and pKa calculations. The first project was developed to overcome this technical limitation: we used an available PB solver in a new approach to perform pKa calculations taking into account the system periodicity. Secondly, to simulate highly charged membranes, a proper treatment of the ionic strength is crucial and a full neutralization of the system is probably too rough an approximation. Hence, we developed a PB-based method to determine the number of ions that should be added to the simulations. With these two problems solved, we developed a new constant-pH molecular dynamics method to deal with charged lipid bilayers (CpHMD-L). This method allows a proper treatment of the periodicity and ionic strength in model membranes. Finally, we applied our methods to three model systems to illustrate the importance of taking into account protonation / conformation coupling in molecular dynamics simulations, in particular when looking at pH dependent phenomena. To the best of our knowledge, this is the only available method that can deal simultaneously with pH considering the protonation / conformation coupling, periodic boundary conditions in protonation free energy calculations, and a careful treatment of ionic strength. This represents a significant improvement in simulations of model biological membranes. It is now possible to move a step forward in the direction of biological membranes since we are now able to simulate a lipid mixture with several different lipids

Modeling Electrostatics and Geometrical Quantities in Molecular Biophysics Using a Gaussian-Based Model of Atoms

Electrostatic and geometric factors are critical to modeling the interactions and solvation effects of biomolecules in the aqueous environments of biological cells as they respectively influence the polar and non-polar components of the associated free energies. Conventional protocols use a hard-sphere model of atoms to devise and study the underlying thermodynamics. But this traditional model tends to overlook some of the important biophysical aspects at the cost of oversimplification of the representation of the solute-solvent environments. Here an alternative and physically appealing model of atoms – a Gaussian-based model, is presented which replaces the hard-sphere model with a smooth density-based description of atoms. This dissertation explains the derivation of a physically appealing dielectric distribution from the Gaussian schematic to model the electrostatics of biomolecules using the implicit-solvent/Poisson-Boltzmann (PB) formalism. It also demonstrates the advantages of using it for computing geometric properties of a molecule such as its volume and surface area (SA) for estimating non-polar portions of the free energy. While highlighting the qualitative importance of the Gaussian-based model, it offers conceptual proofs towards its validity through computational investigations of explicit solvent simulations. It also reports the key features of the Gaussian-based model, which impart to it the capacity of accurately capturing the crucial biophysical factors that characterize biomolecular properties, namely – the effect of intrinsic conformational flexibility and salt distribution. The non-triviality of these factors and their portrayal through the Gaussian models are meticulously discussed. A major theme of this work is the implementation of the Gaussian model of dielectric distribution and volume/SA estimation into the PB solver package called Delphi. These developments illustrate the manner in which the utility of Delphi has been expanded and its reputation as a popular tool for modeling solvation effects with appreciable time-efficacy and accuracy has been enhanced