Understanding the options of atom mapping and ligandHybrid

pmx
1

Dear users and developers,
How do atom mapping options for ligand RBFE differ from each other as in where I use (–H2Heavy option, whether hydrogen be morphed into a heavy atom) and I have been unable to get forward and backward convergence in my FEP calculations.

Whereas, when I use the default options and do not pass any option explicitly, I see many atoms being filled both in A and B state with DUM in the prefix which however seem to be virtual sites/dummy atoms. And in the latter case, my simulations seem to have pretty well backward and forward convergence.

And from my understanding with topologies used in amino acid FEP is that, when I am converting a Halogen atom to a Hydrogen, the hybrid topology should have been pretty straight forward rather than using such dummy atoms/virtual sites. However in such case, the aforementioned problem of forward and backward convergence.

The other query I had about is scaling of dummy masses(–scDUMm), dummy angle parameters( --scDUMa) or dihedral parameters(–scDUMd). What is it that this option essentially does?

It would be kind if anyone could explain me the nuances behind it.
P.S: I have been using BAR and equilibrium free energy calculations.
Warm regards,
Pallav Sengupta,
Guwahati, Assam, India

2

H2Heavy will allow mapping hydrogen atoms in one ligand to heavy atoms in another ligand. I recommend not to use this option.

Many dummies are added when the two ligands share little similarity, yet it does not necessarily mean that the convergence will suffer. As you noticed, the convergence is better with the default mapping options.

These options scale down dummy masses, angle and dihedral force constants between real and dummy atoms. This might be required to ensure simulation stability and to prevent dummy atoms from influencing internal geometry of the rest of the molecule.

1 Like
3

Thanks a lot Vytas.
And there is another query I have regarding soft core potentials during transformations,
how does one determine the value of sc_alpha and sc_sigma, and does these values influence efficient transformation and also affect convergence?
Regards,
Pallav.

4

The soft-core parameters are determined from empirical observations to give the smooth alchemical transition and allow for good convergence. There are some standard values that have been benchmarked to behave well. You can look into the Beutler et al, 1994 paper, where several sets of softcore params were evaluated.

1 Like
5

Hi Vytas,
I have been studying RBFE of protein ligand systems.
This is what I have been observing in case of the hybrid ligand with protein where during equilibration and dynamics it affects the geometry of the nearby molecule, what should have been the ideal case with scaling of dummy masses(–scDUMm), dummy angle parameters( --scDUMa) or dihedral parameters(–scDUMd)? Should the numbers be 1.0 for all these options?

And another problem I have been facing lately is that should I generate the topology files for the protein too using pdb2gmx with the mutated force field or should it be fine with the normal force field i.e CHARMM36 and later add the mutated topologies as done in normal protein-ligand simulations? I mean how do you proceed in such cases.

6

What do you mean “it affects the geometry of the nearby molecule”? Need a clearer formulation of this question.

I usually do not scale bonded parameters and masses.

Regarding the force field, you can use any force field you find suitable. Force fields shipped with pmx are there only for amino acid and nucleic acid mutations.

7

I mean the ligand has stacking interaction with another molecule in the binding pocket of the protein, whereby the hybrid molecule affects the geometry of that molecule itself during the equilibration phase.

Okay.

Fine then.

8

The question then is: does this interaction differ between non-hybrid ligand simulation and hybrid ligand simulation?

9

Yes it is.

10

Just to make sure, does the ligand retain that interaction in a vanilla MD simulation?

If so, you can try scaling down dummy dihedral interactions.

11

what does vanilla MD simulation mean?

In normal classical protein ligand simulation the stacking interaction is maintained all throughout the interaction.

should it be 1.0 as default or one has to try and play?

12

By vanilla I meant without the hybrid topology.

By default the interactions are not-scaled, try scaling them down.

1 Like
13

okay. got it.

okay, thanks.

14

should it be an arbitary number or go through iterations of hit and trial. I have not yet understood the scaling down of these dummy atoms interactions. Is there any reference paper for the same?

15

It is not published, you will need to try it out.

1 Like
16

Have you in any way tried using dual topology for RBFE?

17

Hi Vytas,
I tried it and it does helps. Now the question is to what extent one scales the dummy masses or the angles or the dihedrals? I used a pretty large number(100) for masses, angles and dihedrals and everything seems to be fine with it. But then it would be unrealistic for a dummy atom to have a mass of 100.008.

And for the simulations I had encountered LINCS error, I also used https://ask.bioexcel.eu/t/relative-constraint-deviation-after-lincs/2496/3?u=s.pallav all-bonds in constraints and the simulation seems to be fine.
Regards,
Pallav

18

That requires restraints similar to those in absolute free energies. For that case, I would then directly calculate ABFE

1 Like
19

So you chose to up-scale masses, angles and dihedrals. I don’t think this is meaningful. I would follow my earlier suggestion to down-scale dihedral interactions. The current upscaling does not resolve the issue. The fact that you see an expected behaviour, suggests that the true cause might not be related to the bonded interactions of the dummy atoms.

1 Like
20

I apologise. I was unable to understand it. It is clear now. Clearly saw the numbers in the merged.itp file. Thanks a lot.