How to do MDs
Contents
Amber
Prepare the system (Protein(s) only)
Prepare your protein by removing all hydrogens (like -ignh for pdb2gmx) and feed it to tleap
$cat build.leap source leaprc.protein.ff14SB source leaprc.RNA.OL3 source leaprc.DNA.OL15 source leaprc.water.tip3p mol = loadpdb mol.pdb solvatebox mol TIP3PBOX 10 addions mol K+ 154 Cl- 91 saveamberparm mol ionized.parm7 ionized.rst7
Points to be caution and ALWAYS save you some time:
- Properly distinguish all the components using chain, e.g. chain A for protein subunit A, chain B for DNA strand 1, chain C for DNA strand 2 etc.
- Consider a fix.tcl (see below)
- If a PDB file can go through pdb2gmx, it will be probably fine for tleap
Gromacs
Prepare the system (Protein(s) only)
Justin's tutorials are always the good starting points (here), especially the "Lysozyme in water" one (here).
Simple steps to prepare a system, and it usually works! (script can be found here)
# Create a cubic box of 2 nm on each side and place the system to the center and align its principle axes to the reference axes # -princ usually helps to reduce your system size gmx editconf -f conf.gro -o editconf.gro -bt triclinic -d 2 -c -princ # Solvate the box gmx solvate -cp editconf.gro -o solvate.gro -p topol.top # Add ions to make it neutral and of 0.15 M NaCl # If you want KCl, add -pname K gmx grompp -f ions.mdp -c solvate.gro -o ions.tpr -p topol.top gmx genion -s ions.tpr -o ionized.pdb -conc 0.15 -neutral -p topol.top
Things to note:
- You can find ions.mdp here.
- Biological systems might require surprisingly large boxes, if you are interested, read this.
- You can consider using -d 1.2 if you are confident that your system does not extend its dimension a lot.
Standard MD (Protein(s) only)
A standard MD includes steps as follows:
- Energy minimization
- Annealing
- Pre-equilibrium
- Production run
Here are some standard run scripts.
gmx mdrun -pin on -ntmpi 1 -ntomp 12 -v -deffnm
Things to note:
- In step2-3, position restraints have been added to favor the initial structures (usually crystal)
- You may consider release the restraints step-wisely
Trajectory manipulation
Usually you will have multiple trajectories from production runs, you would like to combine them and remove the PBC effects in order to calculate something or make videos.
gmx trjcat -f *.xtc -o xtc -dt 100 -n ndx
Notes:
- -settime if you did not continue
gmx trjconv -f pre.trr -o unwrap.xtc -s tpr -n ndx -pbc mol (-dt 100) gmx trjconv -f unwrap.xtc -o fit.xtc -s tpr -n ndx -fit rot+trans
Notes:
- Usually protein for least squares fit and centering
- Try "-pbc nojump" if "mol" does not work
- "-pbc nojump" is good for dimer; "-pbc mol" is good for multimer
- !!!Never use "-center" with "-pbc mol"!!!
Split
echo 0 | gmx trjconv -f ed.trr -s ../ed.tpr -split 20 for ii in $(seq 0 10); do jj=$(($ii + 1)) echo 0 | gmx trjconv -f trajout$ii.xtc -s ../ed.tpr -skip 20 -o $jj.gro done rm -f trajout$ii.xtc \#*
Anton
Alt 1:
Load dtr in vmd
pbc unwrap -all animate write trr
gmx trjcat -o xtc gmx trjconv -pbc nojump -o unwrap.xtc
Alt 2:
Load dtr in vmd
animate write trr
Load trr in vmd
pbc unwrap -all animate write trr
gmx trjcat -o unwrap.xtc
Alt 3:
Load dtr in vmd
animate write trr
gmx trjcat -f *.trr -o xtc -cat gmx trjconv -f xtc -o xtc -timestep 120 gmx trjconv -s tpr -pbc res -center -o unwrap.xtc