Help with FES interpretation from OPES

ANTONIO LAUS

unread,

Feb 8, 2024, 6:55:06 AM2/8/24

to PLUMED users

Dear Plumed Community,
I am having a lot of difficulty with FES rescaling of this ligand-protein complex.
I have used an OPES with funnel and two collective variables as my sampling technique.
I ask for your help in interpreting the following graphs to see if I am doing something wrong and if so, suggest me some better way.
My main difficulty is in understanding if I am reweighting the FES right and if there is convergence.

I tried the FES_from_reweighting script first:
FES_from_Reweighting.py --kt 2.4943 --sigma 0.002 --blocks 5 --deltaFat 2.5 --min 1 --max 3.3
Receiving the following error:
"RuntimeWarning: divide by zero encountered in log
fes=-kbt*np.log(np.average(np.exp(-block_fes/kbt),axis=0,weights=safe_block_weight))
average FES uncertainty is: inf"

While the following script works:
FES_from_State.py --kt 2.4943387854 --deltaFat 2.5 --min 1,-5.5 --max 3.5,4.6 --all_stored
But I am not sure if I am doing it right, so I am posting the graphs I have generated from my data and asking if it is possible to figure out, somehow, if there are any errors or inconsistencies.

Thanks in advance to everyone!
Antonio

rct.png

st_evo.png

sigma_d1.png

scatter_evo.mp4

neff_nker.png

free_energy_d1.png

scatterplot.png

d1evo.png

deltaF_evolution.png

free_energy_st.png

zed.png

sigma_st.png

opes_bias.png

Michele Invernizzi

unread,

Feb 9, 2024, 10:26:35 AM2/9/24

to plumed...@googlegroups.com

Dear Antonio,

from your d1evo.png one can see that the simulation never goes back to the binding state, thus it cannot be at convergence. you need multiple back and forth transitions for that. Ideally each block in your block average should contain more than one full transition, otherwise block averaging does not make much sense. You can get a 1D FES estimate for your simulation with FES_from_Reweighting.py --kt 2.4943 --sigma 0.1 --deltaFat 1.7 --cv d1 but it will not be reliable at all.

My suggestion is to try to improve your CVs, or run longer / multiple simulations. I would do these tests using OPES_METAD_EXPLORE, as it is more forgiving than OPES_METAD when the CV is not great.

Best,

Michele

--
You received this message because you are subscribed to the Google Groups "PLUMED users" group.
To unsubscribe from this group and stop receiving emails from it, send an email to plumed-users...@googlegroups.com.
To view this discussion on the web visit https://groups.google.com/d/msgid/plumed-users/e4174bfc-05a7-4a5c-a6c5-f68fed25935an%40googlegroups.com.

Antonio Laus

unread,

Feb 9, 2024, 10:52:02 AM2/9/24

to PLUMED users

Dear Michele,

thank you very much for your reply, continuing with the dynamics I realized that it is returning to the binding position, however I actually did not expect such a long simulation, so I will try to improve the CVs.
I would like to ask you another question, why in your suggestion you reccomend me to set the sigma to 0.1, it is something that I am not clear, because from the kernel file and the graph I notice that the final kernels have an average of 0.002.
Last thing and I conclude, I am leaving the sigma adaptive, in your opinion in the plumed.dat file would it be appropriate to set an accurate sigma, if yes why would it help me?

Sorry for the trouble and thanks again,

Best,

Antonio

Michele Invernizzi

unread,

Feb 12, 2024, 3:17:17 AM2/12/24

to plumed...@googlegroups.com

Dear Antonio,

Choosing the proper width for kernel density estimation (KDE) is a well studied problem, you will find a lot of good material online (e.g. https://aakinshin.net/posts/kde-bw/). When I have to do KDE, especially in 1D, I often use https://docs.scipy.org/doc//scipy-1.4.1/reference/generated/scipy.stats.gaussian_kde.html It's only drawback compared to my script is that it does not handle periodic CVs.

The adaptive sigma should be a good default, especially for OPES_METAD_EXPLORE, but for 1D systems it's also not too difficult to choose SIGMA manually, in a similar way as for KDE. It might give you more consistent performance across multiple runs.

In my experience, the biggest bottleneck is the choice of the CVs. Fine tuning the other parameters might give you some speedup, but if you don't see any recrossing with the default OPES_METAD_EXPLORE (testing a range of reasonable BARRIER values) it's hard to get something out of it.

Best,

Michele

To view this discussion on the web visit https://groups.google.com/d/msgid/plumed-users/d0845e93-7e2c-4af5-b1a1-c9547693e0d6n%40googlegroups.com.

ANTONIO LAUS

unread,

Feb 15, 2024, 6:40:00 AM2/15/24

to PLUMED users

Dear Michele,

thank you very much for your reply , I am using both distance and a collective variable based on the deep-TDA neural network model as collective variables. Indeed it may not be ideal because of the high mobility that the protein has, but I have reason to think that the problem could also be due to an artifact that I notice by prolonging the run, i.e. part of the binding site has a helix that is very distant from the other neighboring helices and this leads the molecule to come out on the other side but because there is the funnel it results as a second site because of the attraction by some amino acid and the funnel wall, I attach the image hoping to be as clear as possible. Ideally, it would be ideal to have a funnel with 2 outputs, I think, or llimiting the RMSD of the protein, but these are just guesses.
In any case, I am doing several tests and as soon as I have more results I will share them here, sorry to bother you but you are the only channel I can compare with.

Best,

Antonio

Screenshot from 2024-02-15 12-25-47.png

Screenshot from 2024-02-15 12-39-00.png

Donkiss Zidane

unread,

Mar 12, 2024, 7:30:28 AM3/12/24

to PLUMED users

Reply all

Reply to author

Forward

Help with FES interpretation from OPES_METAD

ANTONIO LAUS

Michele Invernizzi

Antonio Laus

Michele Invernizzi

ANTONIO LAUS

Donkiss Zidane