# Umbrella Samplig Protein Ligand Gromacs

**URL:** <https://gromacs.bioexcel.eu/t/umbrella-samplig-protein-ligand-gromacs/6094>\
**Category:** User discussions\
**Created:** [April 3, 2023, 12:59pm UTC](https://gromacs.bioexcel.eu/t/umbrella-samplig-protein-ligand-gromacs/6094 "2023-04-03T12:59:25Z")\
**Posts on this page:** 7\
**Page:** 1

<div class="post-metadata">

**Author:** ![alex789](https://avatars.discourse-cdn.com/v4/letter/a/bc8723/32.png) [@alex789](https://gromacs.bioexcel.eu/u/alex789)\
**Post date:** [April 3, 2023, 12:59pm UTC](https://gromacs.bioexcel.eu/t/umbrella-samplig-protein-ligand-gromacs/6094/1 "2023-04-03T12:59:25Z")

</div>

GROMACS version:2021  
GROMACS modification: No

Hello,I am used in performing MD protein-ligand simulations in Gromacs and now I would like to learn umbrella samplig. I have tried this tutorial [Umbrella Sampling](http://www.mdtutorials.com/gmx/umbrella/06_umbrella.html)  
but it isn’t about protein-ligand systems, so I don’t have proper results. Has anyone any recommendations or examples about this topic?

---

<div class="post-metadata">

**Author:** ![Karis](https://avatars.discourse-cdn.com/v4/letter/k/ed8c4c/32.png) [@Karis](https://gromacs.bioexcel.eu/u/Karis)\
**Post date:** [April 3, 2023, 2:43pm UTC](https://gromacs.bioexcel.eu/t/umbrella-samplig-protein-ligand-gromacs/6094/2 "2023-04-03T14:43:28Z")

</div>

Hi Alex, as a starting point, if you were able to do the tutorial all the way through and your results match the tutorial’s graphs in the data analysis section, that should be a pretty good indicator that you’re on the right track. Is there a particular step that you’re unclear on or is causing issues with the simulation?

---

<div class="post-metadata">

**Author:** ![alex789](https://avatars.discourse-cdn.com/v4/letter/a/bc8723/32.png) [@alex789](https://gromacs.bioexcel.eu/u/alex789)\
**Post date:** [April 3, 2023, 4:04pm UTC](https://gromacs.bioexcel.eu/t/umbrella-samplig-protein-ligand-gromacs/6094/3 "2023-04-03T16:04:34Z")

</div>

I have done the tutorial firstly.  
My problem is in the pulling and then in the configurations.  
I tried pulling using the same options as in the tutorial, but as group 1 the protein and as group the ligand and also the pulling to be done on the y axis (I have created a big box, bigger on y axis, similar to tutorial). However, in the xtc file, the ligand moved insignificantly, definetely not 5 nm. On the configuratins the [http://www.mdtutorials.com/gmx/umbrella/Files/get\_distances.sh](http://www.mdtutorials.com/gmx/umbrella/Files/get_distances.sh) script couldn’t calcuate the distances, my problem is that in the tutorial it is not named which is the index.

Thank you for your time

---

<div class="post-metadata">

**Author:** ![Karis](https://avatars.discourse-cdn.com/v4/letter/k/ed8c4c/32.png) [@Karis](https://gromacs.bioexcel.eu/u/Karis)\
**Post date:** [April 3, 2023, 6:24pm UTC](https://gromacs.bioexcel.eu/t/umbrella-samplig-protein-ligand-gromacs/6094/4 "2023-04-03T18:24:05Z")

</div>

To clarify, you’re having difficulty pulling the ligand from the protein by the desired distance, and can’t get the get\_distances.sh script to work due to the script using different indices, is that correct?

On pulling the ligand, different proteins might need different forces to pull apart. The simulation pulls the protein by attaching a spring, then pulling on that spring at a certain rate. If you want to edit the mdp file to get a larger distance, you can increase the force by either increasing the stiffness of the spring (pull\_coord1\_k) or increasing the rate at which you pull (pull\_coord1\_rate). You could also run the simulation for a longer period of time (nsteps), so the spring has more time to pull the protein and ligand apart. Keep in mind that while increasing the force could pull apart the ligand and protein faster, the protein under high force might not behave as it would in a real physical system, so how you adjust the system is gonna depend on what you want out of the simulation.

On the index used, when you made index.ndx, did you name the ligand and protein something other than ‘Chain\_A’ and ‘Chain\_B’? get\_distances.sh is calculating distances for ‘com of group “Chain\_A” plus com of group “Chain\_B”’, so if you change the names of the groups in index.ndx and get\_distances.sh to match, then it should work.

---

<div class="post-metadata">

**Author:** ![jalemkul](https://dub1.discourse-cdn.com/flex017/user_avatar/gromacs.bioexcel.eu/jalemkul/32/18_2.png) [@jalemkul](https://gromacs.bioexcel.eu/u/jalemkul)\
**Post date:** [April 3, 2023, 6:39pm UTC](https://gromacs.bioexcel.eu/t/umbrella-samplig-protein-ligand-gromacs/6094/5 "2023-04-03T18:39:12Z")

</div>

I can quickly weigh in that a protein-ligand complex is quite different from the tutorial.

One should define the pulling vector (reaction coordinate) as some group of residues in the binding/active site and the ligand itself, such that the vector connecting them projects out along a path that exits the binding/active site. This is non-trivial and highly dependent upon the system. It may also not be possible for all systems.

The bias should be applied in all 3 spatial dimensions, not just along one as in the tutorial (which is a very specific case in which this works).

---

<div class="post-metadata">

**Author:** ![alex789](https://avatars.discourse-cdn.com/v4/letter/a/bc8723/32.png) [@alex789](https://gromacs.bioexcel.eu/u/alex789)\
**Post date:** [April 3, 2023, 7:58pm UTC](https://gromacs.bioexcel.eu/t/umbrella-samplig-protein-ligand-gromacs/6094/6 "2023-04-03T19:58:40Z")

</div>

Thank you for help.  
Forgot to mension, my system is from molecular docking.  
If i get it correct, I should define on chain A, a residue that is on the way of the ligand to exit from the binding site and on chain B the ligand. Also, I should change the distance into direction in pull\_coord1\_geometry. Apart from that, should I make any other change in the pull code?

; Pull code  
pull = yes  
pull\_ncoords = 1 ; only one reaction coordinate  
pull\_ngroups = 2 ; two groups defining one reaction coordinate  
pull\_group1\_name = Chain\_A  
pull\_group2\_name = Chain\_B  
pull\_coord1\_type = umbrella ; harmonic potential  
pull\_coord1\_geometry = direction;  
pull\_coord1\_dim = Y Y Y  
pull\_coord1\_groups = 1 2  
pull\_coord1\_start = yes ; define initial COM distance \> 0  
pull\_coord1\_rate = 0.01 ; 0.01 nm per ps = 10 nm per ns  
pull\_coord1\_k = 1000 ; kJ mol^-1 nm^-2

---

<div class="post-metadata">

**Author:** ![Karis](https://avatars.discourse-cdn.com/v4/letter/k/ed8c4c/32.png) [@Karis](https://gromacs.bioexcel.eu/u/Karis)\
**Post date:** [April 5, 2023, 4:21am UTC](https://gromacs.bioexcel.eu/t/umbrella-samplig-protein-ligand-gromacs/6094/7 "2023-04-05T04:21:29Z")

</div>

To my understanding, if you’re using pull\_coord1\_geometry = direction, you also need to define pull\_coord1\_vec ([https://manual.gromacs.org/documentation/2018/user-guide/mdp-options.html#mdp-pull-coord1-vec](https://manual.gromacs.org/documentation/2018/user-guide/mdp-options.html#mdp-pull-coord1-vec)) such that the pulling is in the direction the ligand takes to exit the binding/active site.
