Pulling MD - pullforce has multiplie peaks

User discussions
, ,
1

GROMACS version: 2022.3
GROMACS modification: No
Dear all,

I am trying to follow the steps of the umbrella sampling tutorial as provided by Justin Lemkul for a separate system. Following the tutorial with the initial protein provided in the tutorial works as expected so I expect some issues with my specific system:

I have simulated the adsorption of a peptide to a polymeric surface (extending in x-z direction) using the AMBER99SB-ILDN forcefield. I have extracted the final frame of the adsorption simulation trajectory and can show that peptide and polymer are in close proximity to each other (multiple residues of the peptide are within 3 A of the polymer surface). Using this final frame I wanted to calculate the free energy of binding using the pulling / umbrella sampling approach as described in the tutorial above.
I placed the adsorbed peptide-polymer structure in a box with sufficiently increased y-coordinate (height of peptide-polymer complex: 4.8 nm, pulling distance: 5 nm; y of new box: 15 nm) and solvated the enlarged box and added ions. Naturally, I chose y as the reaction coordinated for the pulling simulation.
I minimized and equilibrated the system prior to pulling as described in the tutorial. This is my pull parameter input file:
md_pull.mdp (2.1 KB)
This is the input file for npt before umbrella sampling:
npt_umbrella.mdp (2.0 KB)
And this is my umbrella sampling input parameter file:
md_umbrella.mdp (2.0 KB)

There is no issue during the simulation. The peptide is pulled away from the surface as expected. The issues arise during analysis. First and foremost I notice that my pulling force never seems to form one distinct peak. Instead, it forms multiple peaks and it appears that the detachment of the peptide from the polymer surface does not impose a larger energy barrier than pulling the peptide through solvent water.

image
image969×754 8.67 KB

Out of curiosity, I used this pulling trajectory for umbrella sampling. Even with a sampling window of 0.1 nm some bins are still poorly sampled (e.g. Warning, poor sampling bin 190 (z=0.9095).). Hence, of course the free energy of binding calculation is not performed over the complete reaction coordinate distance but stops at the poorly sampled reaction distance. However, up until that distance of 0.9095 nm the PMF also never exceeds 1 kJ/mol. Just for comparision, MMPBSA analysis reveals a free energy of binding for this complex of roughly -30 kJ/mol.

Am I wrong in my understanding of the free energy of binding / unbinding that PMF and MMPBSA should yield somewhat comparable results for the same complex? Have I setup the pulling wrongly or has anyone have experience with similar simulations and is able to assist me? Are there any changes I can make to

Any help or suggestions are greatly appreciated!

Thank you very much!

2

You don’t write how large your peptide is, but I think that pulling any molecules with a speed of 10 nm/ns would result in a plot similar to what you post. Your pulling speed is at least one and probably more orders of magnitude too high.

I would also think that umbrella sampling is not going to work in practice as it will be extremely difficult to properly sampling positions where the peptide is part of the times adsorbed and part of the time desorbed.

3

Hi, thank you for your response!
The peptide is 770 atoms large. I have now tried 1 nm/ns, 0.5 nm/ns and 0.1 nm/ns pulling speeds, with respectively increased number of steps to cover the same distance. Unfortunately, all pullforce diagrams look similar.

Does it make a difference if the “true” final frame of a adsorption simulation is used (including original water molecules) or if the peptide and polymer are extracted from the adsorption simulation, the simulation box is enlarged and the system is then re-solvated?

4

You define -DPOSRES_B. Does this impose a restraint on any of the molecules you are pulling? Especially Chain_B?

5

After resolvation you should re-equilibrate the system. Apart from that I don’t think there should be any problems.

6

I don’t think this is going to work in practice. I assume the 770 atom peptide is quite flexible. It is then practically impossible to reasonably sample equilibrium states at different distances. The times scales of switching between partly adsorbed and desorbed states will be too long.

7

Thank you for your reply!
-DPOSRES_B is the position restrain for the polymer to which the peptide is adsorbed. I want to pull the peptide away from the polymer hence I supposed that this would be the way to do it.
After resolvation, I re-equilibrated the system, so that should not be the problem.

8

Thanks for the reply! During the pulling simulation, the peptide actually remains quite stable in its adsorped confirmation.

I found another approach in literature (https://www.pnas.org/doi/10.1073/pnas.0707879105) where the spring is not attached to the COM of the peptide but rather to an exposed terminus of the peptide. I will try to re-create that approach and see if I can get reasonable results.

Thank you all for your help so far!

9

To follow up on my previous reply, if it is chain B you want as reference you should try setting

pull_group1_name        = Chain_B
pull_group2_name        = Chain_A

And I don’t think you’d need any restraints on chain B.

I’m not saying this will solve potential sampling issues, but I think that your pulling force would look better.