GROMACS version: 2019.3
GROMACS modification: No
Here post your question
Hi there,
I just realised I didn’t position my protein in the centre of the box at the start so the protein came out of the box and reappeared on the other side, which could be indicated from the rmsd result of my MD run (figure shown below), as told by the postdoc in my group.
I was thinking if I could fix this problem without restarting the whole MD again. I tried to use gmx trjconv -s md_0_1.tpr -f md_0_1.xtc -o md_0_1_fit-rot_trans.xtc -ur compact -fit rot+trans to centre the protein and then generate the trajectory files (rmsd.xvg etc.) but didn’t succeed. It just changed its shape:
I was just wondering do you have any suggestions on how to deal with it?
Many thanks.
This doesn’t do any centering. If you only have a single protein (not a complex), you can clean everything up (centering and making molecules whole) with gmx trjconv -center -pbc mol.
Thank you very much for the reply.
However, sadly it doesn’t work. Maybe I’ve used it in the wrong way, but the RMSD result seemed to have even more spikes. These are the codes I used to centre my protein:
gmx trjconv -f md_0_1.xtc -s md_0_1.tpr -center -pbc mol
Firstly select ‘protein’ and then ‘system’
Note that major changes are planned in future for trjconv, to improve usability and utility.
Will write xtc: Compressed trajectory (portable xdr format): xtc
Reading file md_0_1.tpr, VERSION 2019.3 (single precision)
Reading file md_0_1.tpr, VERSION 2019.3 (single precision)
Select group for centering
Group 0 ( System) has 130119 elements
Group 1 ( Protein) has 6600 elements
Group 2 ( Protein-H) has 3336 elements
Group 3 ( C-alpha) has 442 elements
Group 4 ( Backbone) has 1326 elements
Group 5 ( MainChain) has 1770 elements
Group 6 ( MainChain+Cb) has 2179 elements
Group 7 ( MainChain+H) has 2192 elements
Group 8 ( SideChain) has 4408 elements
Group 9 ( SideChain-H) has 1566 elements
Group 10 ( Prot-Masses) has 6600 elements
Group 11 ( non-Protein) has 123519 elements
Group 12 ( Water) has 122700 elements
Group 13 ( SOL) has 122700 elements
Group 14 ( non-Water) has 7419 elements
Group 15 ( Ion) has 819 elements
Group 16 ( NA) has 406 elements
Group 17 ( CL) has 413 elements
Group 18 ( Water_and_ions) has 123519 elements
Select a group: 1
Selected 1: 'Protein'
Select group for output
Group 0 ( System) has 130119 elements
Group 1 ( Protein) has 6600 elements
Group 2 ( Protein-H) has 3336 elements
Group 3 ( C-alpha) has 442 elements
Group 4 ( Backbone) has 1326 elements
Group 5 ( MainChain) has 1770 elements
Group 6 ( MainChain+Cb) has 2179 elements
Group 7 ( MainChain+H) has 2192 elements
Group 8 ( SideChain) has 4408 elements
Group 9 ( SideChain-H) has 1566 elements
Group 10 ( Prot-Masses) has 6600 elements
Group 11 ( non-Protein) has 123519 elements
Group 12 ( Water) has 122700 elements
Group 13 ( SOL) has 122700 elements
Group 14 ( non-Water) has 7419 elements
Group 15 ( Ion) has 819 elements
Group 16 ( NA) has 406 elements
Group 17 ( CL) has 413 elements
Group 18 ( Water_and_ions) has 123519 elements
Select a group: 0
Selected 0: 'System'
Reading frame 0 time 0.000
Precision of md_0_1.xtc is 0.001 (nm)
Using output precision of 0.001 (nm)
Back Off! I just backed up trajout.xtc to ./#trajout.xtc.1#
Last frame 10000 time 100000.000 -> frame 10000 time 100000.000
Then I generated my rmsd file with this:
gmx rms -s md_0_1.tpr -f md_0_1.xtc -o rmsd.xvg -tu ns
Selected “1: protein” in the next two steps.
I don’t know in which step I made some mistakes and just wanna post here in case you have any thoughts. Hope the questions all make sense.
Here, you did not specify an output file with -o so you will get the default trajout.xtc.
Here, you analyzed the old trajectory that you did not center. You should perform the analysis on trajout.xtc (or, ideally, something with a more descriptive file name).
Hi,
Thank you very much for the reply.
I followed what you suggested these days but sadly it still didn’t work:
gmx rms -s md_0_1.tpr -f trajout.xtc -o rmsd_trajout.xvg -tu ns
I still got spikes in the RMSD results:
As indicated in the figure, the spikes are different from the previous figure I posted several days ago. My supervisor suspected it was generated because the protein wasn’t structurally aligned with the before-simulation protein, so I just tried to structurally align my protein to the reference one. These are the commands I’ve used:
gmx trjconv -s md_0_1.tpr -f md_0_1.xtc -o md_0_1_noPBC.xtc -pbc mol -center (selecting 1 - Protein and 0 - system to re-center)
gmx trjconv -s md_0_1.tpr -f md_0_1_noPBC.xtc -o md_0_1_align.xtc -fit rot+trans (selecting 1 - Protein and 1 - Protein for output)
gmx rms -s md_0_1.tpr -f md_0_1_align.xtc -o rmsd_align.xvg -tu ns (selecting 1 - Protein and 1 - Protein for output)
I’ve also checked the similar page (Structural Realignment) and tried the commands over there but there’re still spikes in my RMSD results as shown above. Now I started to suspect whether the plotting tool I used for generating RMSD was not appropriate. For your information, I used grace (Grace Home) to plot my RMSD from the rmsd.xvg file. The postdocs in my group also suspected perhaps the rescheduling of the HPC cluster I was using caused this result, but it’s still the same after I did some replicates.
Any ideas would be greatly appreciated, since I’ve basically tried every possible solution I could get from my group as well as the forum.
Many thanks.
Hi. This problem has already been solved and I would like to post the solution here in case anyone else has the similar problem in the future.
It turned out that I didn’t use ‘-merge’ function when I used pdb2gmx during the protein preparation process. My protein is a Fab domain containing two chains and with -merge function will the two chains be merged into one moleculetype, only in which way will not the chains split away from each other during the simulation. Then after the simulation, using the trjconv nojump and rot + trans options to remove translation and rotation of the protein. All is done.
Hi, Mario here. I am new to gromacs, and I have run into such issues, i. e. a dimer that is split during visualization in pymol or vmd. I am wondering if “-pbc nojump” might be enough to visualize the dimer at the centre of the box. I was wondering if “-merge” can alter the md simulation if multiple chains turn into one molecule.
Thanks if you ever reply
The typical approach is to create an index group for one chain and use trjconv to center on it, which will wrap the other chain to it. Simply applying -pbc nojump will not fix this case, but -pbc mol -ur compact -center usually does if you use this approach.
1 Like
sorry to annoy you, but how can one create an index file with a single chain that is not indicated within the .gro files but only in the initial pdb file? I tries to use the option r 1-133 for a single chain, but the number of atoms was higher than that expected
you were absolutely right, I made a custom index, centered it and it worked perfectly and solved the problem for all my simulations, I guess working with indexes is really the way to go with rmsds, it is great we have a forum like this
