GROMACS version: 2020.2
Hi everyone, I’m performing a steered MD run of my ligand through a bilayer. I would like the ligand to move in the z direction only. However, I’m wondering if anyone might know why the ligand drifts (moves in the xy direction in addition to the desired z direction) instead. It is bizarre, since I set pull-coord1-dim to be N N Y, which tells the program that the ligand can only move in the z direction along the vector connecting the two groups.
The commands that I used are shown below,
pull = yes
pull_ncoords = 1 ; only one reaction coordinate
pull_ngroups = 2 ; two groups defining one reaction coordinate
pull_group2_name = BILAYER
pull_group1_name = LIGAND
pull-group2-pbcatom = 3347
pull-group1-pbcatom = 34099
pull-pbc-ref-prev-step-com = yes
pull_coord1_type = umbrella ; harmonic potential
pull_coord1_geometry = distance ; simple distance increase
pull-coord1-dim = N N Y
pull_coord1_groups = 1 2
pull_coord1_start = yes ; define initial COM distance > 0
pull-coord1-init = 0
pull_coord1_rate = -0.0001 ; 0.0001 nm per ps = 0.1 nm per ns
pull_coord1_k = 1000 ; kJ mol^-1 nm^-2
Hi,
This set-up will indeed pull the ligand along the Z-direction only, but that does not mean that your ligand can no longer move freely in the XY-plane. If you want to restrain the XY-position of your ligand, you can add a second pullcoord with dim Y Y N, and a zero pull-rate. However, I don’t know what the advantage would be of fixing the XY-position.
Kind regards,
Wouter
Hi Wouter,
Thanks! Correct me if I am wrong, but is it true that the bilayer shouldn’t be moving, since I would want a straight path along the Z-direction for the ligand to pass through the bilayer from the water phase?
Also, I am a little confused. If I set the bilayer to be the second pullgroup, with dim Y Y N, and a zero pull-rate, would that restrain the XY-position of the bilayer? Could I simply set dim for the bilayer to be N N N?
Thanks so much!
Hi,
I made a typo in my previous message, I meant to say that you could create a second pullcoord (and not a second pullgroup).
The first pullcoord would be identical to the one you have, and takes care of pulling in the z-direction. The second pullcoord would fix the XY-position of the ligand COM to the XY-position of the bilayer COM.
So setting:
pull_ncoords = 2
and adding
pull_coord2_type = umbrella ; harmonic potential
pull_coord2_geometry = distance ; simple distance increase
pull-coord2-dim = Y Y N
pull_coord2_groups = 1 2
pull_coord2_start = yes ; define initial COM distance > 0
pull-coord2-init = 0
pull_coord2_rate = 0 ; restrain in place
pull_coord2_k = 1000 ; kJ mol^-1 nm^-2
The bilayer can still move, but with this set-up, the XY-position of the ligand COM will follow the XY-position of the bilayer COM. Setting dim to be N N N will either result in errors, or perform no pulling at all.
If you want a straight passing of the ligand through the membrane solely for visual purposes, doing this in post-analysis would be better. The more positional restraints you introduce, the more ‘unphysical’ your trajectory will become.
Kind regards,
Wouter
Thanks wouter. I will try this out. My eventual aim is to generate the pmf profiles of the ligand permeating the bilayer from water phase. So is it okay for the ligand to drift in the xy direction as well?
Hi,
For the PMF profile along the bilayer normal (z-direction), you would normally not restrain the xy-coordinates, as you could block certain conformations that could contribute to the ensemble average. I
assume that your ligand only drifts in the xy-plane in the water phase, and not (or noticably less) after it enters the bilayer?
You will probably use umbrella simulations, for which there is a nice tutorial available on the gromacs-tutorials website: GROMACS Tutorials
I don’t know your system, and have no experience with such large ligands. If you are not sure, it’s best to check literature for inspiration.
Kind regards,
Wouter
1 Like