During PBC Removal for my protein in water simulation, I noticed that the protein itself is rotating a lot. I wanted an aligned trajectory for downstream analysis hence I ran
gmx trjconv -s md1.tpr -f md1_compact.xtc -fit rot+trans -o md1_compact_fit.xtc
I chose backbone for fitting. The rotational motion of the protein is no longer present, but the entire water box is now spinning around like crazy and parts of it leaves the pbc box as well(visualized on VMD)!!
Has someone faced this issue before? My problems are compounded by the fact that I need to analyse the solvent properties as well
Dear @abelxf ,
Nothing is really leaving the simulation box (see PBC). The rot+trans fit is basically removing rotational and translational motion of the group you specify, which means that the protein will be nearly immobile at the centre of the box, but this is only a change in the frame of reference you are using. If you think about it, it is like you were sitting on the protein and taking a look around you. From the frame of reference of the box, the box orientation is fixed and you see the protein roto-translate. If now you sit on the protein, you are standing still (the protein is moving with you) and it’s everything else that moves around. That’s what you are looking at, the protein is ‘standing still’ and the rest of the box is moving around it.
Regarding the calculation of solution properties, in principle changing the representation of your coordinates doesn’t change the content of them. In practice, you may want to pay particular attention to how the quantities are calculated, e.g., water diffusivity can be extremely wrong if have PBC jumps which will skew your MSD calculations.
Dear @obZehn,
Thanks a lot! What you have described sounds extremely convincing, but in a separate trajectory, I have observed the water moving around so much that the protein itself has left the solvent with parts of the solvent being outside the solvent, but the protein is still centered within the pbc box, assuming the pbc box that vmd provides is accurate(I have uploaded an image of the same). Hence I’m a bit sceptical about running water dynamics analysis on this trajectory as the water seems to be spinning around for no good reason
I don’t see really any problem anyway, because you have periodic boundary conditions. Again, this rot+trans fit is a postprocessing of the trajectory. Check the trajectory before the fit, you will see that the protein doesn’t leave the water phase.
You might be able to partially solve the visualization problem by running a first time gmx trjconv and specifying the flags -pbc mol -center and centre the protein in the box. Then, you run the fit rot+trans on the output trajectory of trjconv.
Regarding the analysis of water properties, in principle I would run the analysis on the original trajectory. Per se, changing the frame of reference shouldn’t change the properties you measure, but again, it depend on how the algorithms are implemented (e.g. MSD can be horribly wrong if you have jumps in the box, same as RMSD, etc.). Take a look at what properties you want to calculate, the tool to calculate them, and their eventual pitfalls.
Hope this helps!
Dear @obZehn
The steps I had followed for pbc removal prior to use of -fit rot+trans are
- gmx trjconv -pbc whole (The protein was broken across the box)
- gmx trjconv -pbc mol -center (The protein was at the edge of the box after -pbc whole)
- gmx trjconv -fit rot+trans (The protein was rotating around a lot, hence messing up downstream analysis)
Before applying fit rot+trans, no part of my protein was leaving the box. It was safely confined to the center, but rotating a whole bunch. Since for the algorithm I am using to find the group of water molecules of my interest require the protein to be aligned to the z axis(in this case, to the .tpr file), I had to align the trajectory.
Does that seem fine? Thanks a lot for your time
Dear @abelxf
To me this looks fine, as long as it is a matter of visualization. Then, for analysis, it depend on the tool you are using and what the tool expects (e.g. to have unwrapped pbc coordinates etc). If you need the protein aligned along the z axis then I would say that you may also try the gmx editconf tool that has some -align flag!
Dear @obZehn
Thanks a lot! I am actually using an in house python script to identify water within the cavity, which involves alignments, hence my concern regarding the solvent properties. Thanks a lot for your comments!