[gmx-users] about periodic image minimum distance

Justin A. Lemkul jalemkul at vt.edu
Tue Sep 20 13:57:23 CEST 2011



Anna Marabotti wrote:
> Dear gmx-users,
> I'm doing several simulations of a monomeric protein with different 
> ligands. I had several frustrating experiences because I found that at 
> the end of simulations, I found that my complexes violate the minimum 
> periodic image distance.
>  
> I'm using Gromacs 4.5.3 and for all my simulations I set a triclinic box 
> with the command:
> editconf -f prot.gro -o prot_box.gro -bt triclinic -d 1.5 -c
> I think that this dimension is sufficient since I'm using rcoul= 1.0 and 
> rvdw=1.4 nm. I'm currently using amber99sb force field with tip3p water; 
> the system was neutralized with Na.
> I have minimized (to Emtol=500) my system, equilibrated with 20 psec 
> NVT+100 ps NPT (I checked the energies, and they were all constant), 
> then submitted to a full MD for 30 ns.
>  
> When I had a look at the trajectory using VMD, I saw that the protein 
> goes partially out of the box, therefore I used:
> trjconv -f prot_full.xtc -s prot_full.tpr -o prot_fullC.xtc -pbc nojump 
> (selection: System=0)
> and I analysed the new trajectory with:
> g_mindist -f prot_fullC.xtc -s prot_full.tpr -od prot_fullC_mindist.xvg 
> -pi (selection: Protein=1)
>  
> The graph shows that the minimum periodic image is < 1.4 nm for a large 
> part of the trajectories (and even < 1.0 nm in some cases), therefore 
> I'm assuming that the protein interacts with its periodic images (and 
> that the trajectories must be thrown away...)
>  
> I would like to ask some questions. Maybe they are trivial, but I hope 
> the answers will be useful not for me only.
> 1) Am I doing correctly the analyses? Am I indicating correctly the 
> System (option 0) to perform the nojump option, and the Protein (option 
> 1) to analyse the minimum distance?

I would try running trjconv -center on prot_fullC.xtc, centering on the protein, 
to guarantee that you have a fixed reference point.  In reality, this is a 
non-issue; there is no "center" of a periodic system, but it is far easier for 
us to visualize and conceptualize.  That way, if you have a frame that has a 
small minimum image distance, you can easily see what's going on in the context 
of the protein's location within the box.

> 2) What could be the causes for these recurrent periodic image 
> violations in my system(s)? Why does a protein "go out" of the box? Can 
> I argue that the box is too small? Can I argue that the system is not 
> fully relaxed before the production MD?  What are the analyses that I 
> can make to assess these conditions before doing the MD?

There is no "outside" or "inside" in a periodic system.  All molecules can 
traverse boundaries and come back in from the opposite periodic face.  This is 
not a problem.

> 3) What are the solutions? Do I have to increase the box dimensions? Do 
> I have to change the box type? In another case, I solved the problem 
> using a rhombic dodecahedron box, so is this a problem recurrent only 
> for triclinic box, or not? Can a longer stabilization of the system 
> avoid this problem?
>  

If you have a globular protein, you should always use a pseudo-spherical unit 
cell, since their symmetries are complementary (and the simulation will be faster).

I can see no reason a priori why a box-solute distance of 1.5 nm would be 
insufficient, provided that your protein stays stably folded.  I've seen so many 
complaints about g_mindist at this point that I'm beginning to wonder if there's 
a problem with the program rather than what the end user is doing.

> In summary, I'd like to know if there is a way to predict the risk of 
> having a similar problem looking at the system BEFORE simulation is made 
> (and time is wasted...), and if it depends on some incorrect setting of 
> the system (maybe some .mpd option?), or not.
>  

Provided that the protein stays folded, a box-solute distance of no less than 
half the longest cutoff should usually avoid this problem (depending on the 
stability of the system under pressure coupling, so r/2 is an absolute minimum 
and should not actually be used).  Your use of -d 1.5 should preclude any 
issues, although I realize the output of g_mindist contracts that assertion.

-Justin

> Thank you very much for your answers.
> Anna
>  
> ____________________________________________________
> Anna Marabotti, Ph.D.
> Laboratory of Bioinformatics and Computational Biology
> Institute of Food Science, CNR
> Via Roma, 64
> 83100 Avellino (Italy)
> Phone: +39 0825 299651
> Fax: +39 0825 781585
> Email: anna.marabotti at isa.cnr.it <mailto:anna.marabotti at isa.cnr.it>
> Skype account: annam1972
> Web page: http://bioinformatica.isa.cnr.it/anna/anna.htm
>  
> "When a man with a gun meets a man with a pen, the man with a gun is a 
> dead man"
> (Roberto Benigni, about Roberto Saviano)
>  
> 

-- 
========================================

Justin A. Lemkul
Ph.D. Candidate
ICTAS Doctoral Scholar
MILES-IGERT Trainee
Department of Biochemistry
Virginia Tech
Blacksburg, VA
jalemkul[at]vt.edu | (540) 231-9080
http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin

========================================



More information about the gromacs.org_gmx-users mailing list