Tutorial: Simulating AdK with Gromacs
In this tutorial we will perform a molecular dynamics (MD) simulation of the the enzyme adenylate kinase (AdK) in its open conformation and without a ligand bound. We will simulate it in a realistic environment (100 mM NaCl solution at T = 300 K and P = 1 bar) and analyze its structural properties.
For this tutorial we use Gromacs (version 4.5.5) to set up the system, run the simulation, and perform analysis. The overall workflow moves through the following steps:
- Download the tutorial files and organize the work space
- Setup
- obtain structure 4AKE from PDB, select chain A
- use default protonation states
- generate topology
- solvate in water in simulation cell (rhombic dodecahedron)
- add NaCl ions to neutralise and final physiological concentration of 100 mM
- Energy minimization (EM)
- Position restraint equilibration of solvent (MD); NPT (use weak coupling (Berendsen) schemes)
- Equilibrium MD simulation (unrestrained, NPT, use Nosé-Hoover and Parrinello-Rahman temperature and pressure coupling)
- Analysis
- RMSD
- RMSF
- opening (NMP-LID distance)
- radius of gyration
- secondary structure
All input files are provided in the tar file AdKTutorial.tar.bz2 . Start by downloading and uncompressing the file:
curl -O http://becksteinlab.physics.asu.edu/pages/courses/2013/SimBioNano/13/AdKTutorial.tar.bz2 tar -jxvf AdKTutorial.tar.bz2 cd AdKTutorial
Then begin the Simulating AdK with Gromacs Tutorial .
Help with PME-switch in gromacs
My name is Rakesh and I am a PhD student at Molecular Biophysics Unit, Indian Institute of Science. I am following your Adk tutorial on gromacs but I have a doubt regarding the option PME-Switch. From the gromacs mdp options, I see that it is used for constant energy simulations but I also see in forums that this option need not be used especially with the Verlet cutoff in latest versions. So my question is should I use PME or PME-Switch ? If PME-Switch what is the advantage especially with respect to CHARMM27 force field, since I see your mdp file has been optimized for this force field. Also should rcoulomb-switch be specified as I dont see it in your mdp file. What is the rcoulomb and rcoulomb-switch cutoffs I need to use ? Finally you mention in your tutorial for simulation box “for simulations you want to publish this number should be 1.2…1.5 nm so that the electrostatic interactions between copies of the protein across periodic boundaries are sufficiently screened.” I see in some publications that the 0.9nm to 1.0nm being used for CHARMM27 force field in gromacs. Taking into account minimum image convention should it not be greater than 1.2nm or is it highly dependent on the protein you are simulating or rcouloumb cutoff. I would be grateful if you could clear my doubts on these aspects or point me to the relevant literature.The AdK tutorial might need some updating and in particular there’s no guarantee that any of these settings correspond to what is needed for publishable results. The main purpose of the tutorial is to quickly introduce students to running a simple MD simulation with Gromacs. In general you should look in published papers for what settings are acceptable and of course for Gromacs the Gromacs mailing list is an important resource.
For the non-bonded parameters with CHARMM36 and Gromacs 4.6.x we are currently using the Verlet scheme (
cutoff-scheme = Verlet) with PMEcoulomb-modifier = potential-shift-verletandvdw-type = cut-offwithvdw-modifier = potential-shift-verletandvdw = 1.2and the CHARMM TIP3PS model. These settings were chose as those that best reproduced the behavior of POPC and POPE membranes, compared to simulations in CHARMM with native CHARMM36 non-bonded settings. Following the insights from Piggot et al (doi: 10.1021/ct3003157) we also used LINCS with order 4 and 2 iterations. In general, it is not easy to directly transfer settings between different MD codes.However, note that Justin Lemkul recommends a slightly different set of parameters for CHARMM36 and Gromacs 5.x and if in doubt, I would follow his advice.
Given typical salt concentrations of 100 mM and typical Debye lengths of 1 nm we currently consider about 1.2 nm distance between periodic copies of a protein the acceptable minimum but you need to decide for yourself (and possibly convince referees that smaller numbers are appropriate).