I have completed a protein-ligand simulation for 150 ns using Charmm36 force field in accordance with the tutorials written by Dr. Lemkul. Now that my simulation is over, I am interested in checking whether my protein has seen its periodic image during the course of the simulation. I have used a cubic box, (system size: 6.055 7.241 5.252 (nm), diameter: 7.488 (nm), center: 8.346 5.001 4.642 (nm)) with a distance of 1.0 nm (-d 1.0) from the box edge. My vdw cut-off is set to 1.2 nm (rvdw =1.2) and long-range electrostatic cut-off is 1.2 nm (rcoulomb=1.2). Now, I am using the gmx mindist module of gromacs to check the minimum distance between the protein and its periodic image during the course of the simulation. To begin with, I have corrected my trajectory for pbc by recentering the protein and rewraping the molecules within the unit cell by invoking the trjconv command:
gmx trjconv -s md_0_150.tpr -f md_0_150.xtc -o md_0_150_center.xtc -center -pbc mol -ur compact
Next, I used the command:
gmx mindist -s md_0_150.tpr -f md_0_150_center.xtc -n index.ndx -pi -od mindist.xvg
Then selected C-alfa when prompted
Once the calculation got over, the following output was printed on the terminal:
The shortest periodic distance is 2.15047 (nm) at time 144990 (ps), between atoms 109 and 3663.
I am attaching my graph for reference,
My question is am I doing this correctly and how can I interpret the results of gmx mindist according to my system. Is there a chance my protein can interact with its periodic image? I am a beginner in gromacs and need some expert advice in this regard.
Thanks in advance.
Hi,
The distance between the protein and its periodic image should not be smaller than the cut-off used to describe non-bonded interactions. To be sure that this is your case one should considered all proteins atoms (not only the Calpha) when one calculates the minimum distance.
\Alessandra
Thank you very much for your response. I have one small follow-up question. Do I need to perform the minimum distance analysis using the raw trajectory file md.xtc or on the trajectory md_center.xtc file that has been corrected for PBC using the gmx trjconv module with flags -center, -pbc mol, and -ur compact (like I am doing right now). Please give your suggestion. Thanks in advance.
Hi,
it should not matter in the calculation of the minimum distance of one molecule with its image if you use the raw trajectory or you post process it with the option -center -ur -pbc
\Alessandra
My first doubt is related to the output mentioned above. I understand that this is the distance between atoms 131 and 4854. Is it that the atom 131 belongs to the protein from the main simulation cell and atom 4854 is an atom from the periodic image of this protein?
Also, what I understand is that the distance between the protein and its periodic image should not be smaller than the cut-off used to describe non-bonded interactions. Since the cut-off used for bonded and non-bonded simulations is 1.2 nm and the least distance I am getting now is 1.12 nm so do I need to repeat the simulations again? I have used a cubic box, (system size: 6.055 7.241 5.252 (nm), diameter: 7.488 (nm), center: 8.346 5.001 4.642 (nm)) with a distance of 1.0 nm (-d 1.0) from the box edge. What could be the dimensions of my box now if I want to repeat my simulation? Should I also increase the distance of my protein from the edge of the box from -d 1.0 to 1.2?
Hi,
If the protein is globular, I suggest to use a truncated octahedron as a box (see here for details https://manual.gromacs.org/current/reference-manual/algorithms/periodic-boundary-conditions.html?highlight=box or this tutorial Introduction to Molecular Dynamics — GROMACS tutorials https://tutorials.gromacs.org documentation).
I have noticed that your original box was not cubic (a=b=c), and also at the begin of the simulation the minimum distance was close to 1.2. Maybe good to cross check how it was built the original box and if the protein change conformation. By looking at the simulations and by revising how the box will be built, may help you to understand how much to increase -d. e.i. Did your protein structure enlarge during the simulation? Then it may be good to use a -d option on the large side.
I hope it helps
\Alessandra
