GROMACS version: 2021.4
GROMACS modification: No
Cross-posted on Research Gate
I’m running 100 ns long MD simulations of a small peptide in isolation or bound to a partnering protein. My goal is to produce a plot of the residues in double helical form from start to end along with a video of the structure using Chimera.
The steps I take after a MD run are…
Step 1:
Center the trajectory file:
gmx trjconv -s md_0_1.tpr -f md_0_1.xtc -o md_0_1_noPBC.xtc -pbc mol -center
— Type 1 for protein and then on the next prompt 0 for system.
Step 2:
Calculate the back bone rmsd of the entire complex:
gmx rms -s em.tpr -f md_0_1_noPBC.xtc -o rmsd_xtal.xvg -tu ns
— Type 4 then 4 for backbone RMSD calculations
Step 3:
I then make an index file of the two chains:
gmx make_ndx -f md_0_1.gro -o index.ndx
— Type “splitch 1” hit enter, type the number of what chain I want to reference (typically “2”), then type “q” to exit
Step 4:
From there I want to produce a trajectory file (.xtc) from the index and an updated topology (.tpr) file from my MD run using the index file:
Step 5:
Trajectory of just chain 2:
gmx trjconv -f md_0_1.xtc -s md_0_1.tpr -n index.ndx -o chain2_traj.xtc
— Type the chain I want to reference (typically “19”)
Step 6:
Topology of just chain 2:
gmx convert-tpr -s md_0_1.tpr -n index.ndx -o md_chain2.tpr
— Type the chain I want to reference (typically “19”)
Step 7:
I then center the updated files:
gmx trjconv -s md_chain2.tpr -f chain2_traj.xtc -o chain2_traj_noPBC.xtc -pbc mol -center
— Type “1” then “0”
Step 8:
I then perform my analysis using dssp and rms
gmx do_dssp -f chain2_traj_noPBC.xtc -s md_chain2.tpr -sc scount_chain2.xvg -o ss_chain2.xpm -dt 10 -sss H
— Type “5” for mainchain
Step 9:
gmx rms -s md_chain2.tpr -f chain2_traj_noPBC.xtc -o rmsd_xtal_chain2.xvg -tu ns
— Type “4” then “4” for backbone rmsd
*** The issue I run across is when I open the tpr and xtc files in Chimera with the MD viewer, the peptide is completely wrong. What I’ve deduced is that the topology (tpr) file for some reason gives each atom in the first few residues their own residue id’s. Thus residue 1 actually makes up residues 1-9 (or something similar). However, when I produce a pdb from the trajectory file (.xtc) with:
‘’’
gmx trjconv -s md_chain2.tpr -f chain2_traj.xtc -dt 100 -o trj.pdb
‘’’
The peptide is correct. So I believe the issue is that the topology file is being produced improperly due to a bug or something wrong with my steps above. I’ve also tried using the em.tpr file for input on step 6, but I run across the same issue. This is happening with several different peptides. One peptide that is 65 amino acids long and one that is 18 amino acids long. On the 65mer, it appears to only happen on the first 4-5 residues. ***