A Technical Introduction to the MD-OPEP Simulation Tools Jessica Nasica Université de Montréal Montréal, Québec, Canada
Outline • Overview of the simulation method • Why use “simulator”? • The OPEP force-field • Running the simulation: what you need to get started • The simulation step by step • The output of the “simulator” • The analysis tools • The PTWHAM analysis • The RMSD analysis • Secondary structure tools • Contacts tools • Clustering tools • Graphical tools and various useful scripts
The simulation method • REMD simulation ( Replica Exchange Molecular Dynamics ) • A simplified coarse-grained potential: OPEP
Why use “simulator”? • It can be used: • To study aggregation processes • To check for the stability of a particular molecular assembly • The simplified force-field allows to reach bigger time scales more efficiently: • To study long-range proteins motions • To extract accurate thermodynamics properties P. Derreumaux, N. Mousseau – J. Chem. Phys. 126, 025101 (2007)
The OPEP force-field OPEP = Optimized Potential for Efficient structure Prediction • Reduced-protein model: • A 6-particle model with a detailed representation of the backbone (except for Proline) • Each side chain is represented as one particle and defined by one centroid • Implicit solvent solvent effects included in the interaction parameters Maupetit et al. - Proteins 2007; 69:394-408 / Chebaro et al. – J. Phys. Chem. B 2008
The OPEP force-field • OPEP energy function: • Local Forces: Includes changes in bond lengths and valence angles for all particles and changes in improper torsions of the side chains. • Nonbonded Forces: Van der Waals, electric and hydrophobic forces short- and long- range interactions are computed separately. Use of a pairwise contact potential between side chains represented either by a 12-6 potential or by a 6-potential. • Hydrogen-bonding Forces: 2 terms 2-body H-bonds 4-body effects = cooperative energies between H-bonds Maupetit et al. - Proteins 2007; 69:394-408 /Chebaro et al. – J. Phys. Chem. B 2008
Running the “simulator” • What you need to get started: • The executable file after compilation of the code • The parameter file ‘simulator.sh’ • Your ‘.pdb’ file in an OPEP format • The corresponding topology file (.top), residue description file (.list) and ‘ichain.dat’ file • The ‘cutoff.dat’ and ‘scale.dat’ files • A link file and a qsub script
The simulation step by step • Step 1:Minimization: finding a stable starting point It makes sure that the structure’s energy is minimized using ART. At the saddle point: • The 2 phases of ART: • Activation phase the structure is pushed towards a saddle point • Relaxation phase the structure is pushed slightly over the saddle point and relaxed to a new local energy minimum • Mousseau et al. – Frontiers in Bioscience 13 - 2008; 4495-4516
The simulation step by step • Step 1:Minimization: the output The minimized configuration is written in the file “relaxed_conformation.pdb”. The data resulting from each minimization step is written in “log.file”: Total Energy . towards a minimum Net Force to 0 Decrease in potential energy term Velocity Decrease in kinetic energy term Rmsd value Structure is undergoing a configurational change
The simulation step by step • Step 2: Thermalization: It thermalizes the configurations by heating up by stages until it reaches the target temperature. Here, in “simulator.sh”: E = T = -98.5290 0.00 K E1 = -93.859 T1 = 59.25 K E2 = -94.2221 T2 = 88.88 K E3 = -72.7113 T3 = 133.32 K E4 = 25.8615 T4 = 199.98 K E5 = 103.8417 T5 = 299.97 K Initial conformation Relaxedconformation (from minimization step) Thermalization conformations 1 5 Energies are in Kcal/mol
The simulation step by step • Step 3: MD calculation of forces: using the velocity-Verlet algorithm for integrating Newton’s equation of motion for each particle: where, the force on particle i. (Highly simplified description of MD)
The simulation step by step • Step 4: Writing of the .pdb files: At each desired time step, the new configuration is written in the “min” files for each temperature.
Output files The configurations are sorted by temperature.
The PTWHAM tool • Allows a temporal correlation in the data, using autocorrelation analysis, to compute equilibrium averages. we can derive thermodynamical properties, including data from each temperature.
The PTWHAM tool • How to use it: • ./ptwham_first will read data in each pXX/min file, calculate the Rg, rmsd and end-to-end distance, and output files simXX.txt, wham_parameters.dat, beta.dat • edit wham_parameters.dat • ./ptwham_second mintime maxtime writes “averages.txt” containing U (total energy), rmsd, Cv, end-to-end distance, writes “free_energies.txt” containing Uk (average energy per replica), free energy and entropy and writes “deviations.txt” containing all the uncertainties on the above observables.
The PTWHAM tool • The output: Using 2 python scripts, we obtain the following plots:
The RMSD tool • Calculates the rms distance between 2 configurations or between 1 configuration and a list of configurations taking into account their individual clusters for more accuracy. • Recognizes clusters from hydrogen bonds formed using the DSSP definition of a hydrogen bond: the DSSP algorithm identifies a H-bond if E<-0.5 Kcal/mol. • The program searches all the peptides forming H-bonds and rearranges the pdb file accounting for the chain order in the clusters formed. How to use it: rmsd conf1.pdb list It needs “ichain.dat”.
The RMSD tool The output files are : • rmsd.txt containing all the rmsd values for each configuration in the list • id_order.txt containing all the clusters formed for each configuration in the list • parallel_planes.txt containing a list of all the clusters being parallel to each other
The Replicas tool • A python script that generates a graph showing acceptance probability for replica exchange. • Need “log.file” and “replicas.dat”
Additional useful scripts • Pymol scripts allowing to prepare the pdb files and make movies out of the configurations unique_pymol command: unique_pymol file skip • Python scripts: • extract_confs extracts a subset of configurations from a single file containing a list of configurations. • repair_min repairs min files that are broken by a crash, removing all the incorrect lines. • Graphics scripts: • trace_averages.py ( WHAM analysis ) command: ./trace_averages.py “averages.txt” “title” “1” • trace_freeenergies.py ( WHAM analysis ) command: ./trace_freeenergies.py “free_energies.txt” “title” “1”
Where to find our packages? You will find an updated version of the analysis tools on the wiki:http://riel.pmc.umontreal.ca/groups/biophysique/