Funnel & OPES advice and reweighting

[ANTONIO LAUS]

Hi everyone, I would like to try using funnels with opes. I would like to get some advice regarding the coupling of these two techniques both in terms of setup and reweighting.
- As far as the setup is concerned, I would like to understand if it would be enough to set up funnel and simply replace the metadynamic parameters with those of OPES, or are there any special tricks I should pay attention to?
- As far as reweighting is concerned, since the OPES does not output a HILLS file, I would like to understand if it is possible to do this canonically, i.e. as I would do with the funnel metadynamics and then with ffs_tk or the like, or should I use another method?

Thank you to everyone who answers,
Antonio

[Michele Invernizzi]

Hi Antonio,

There is no problem in combining OPES with a funnel-shaped fixed bias potential, and it has been done in the past, see e.g. https://www.nature.com/articles/s41467-020-20310-0 they also provide all input files.

To compute the reweighting (as with plain metadynamics) you do not need the HILLS file. You can find suggestions on how to perform post processing with OPES in this tutorial https://www.plumed.org/doc-v2.9/user-doc/html/masterclass-22-03.html, probably the script that can help you most is this one https://github.com/invemichele/masterclass-22-03/blob/master/scripts/FES_from_Reweighting.py

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/be277de8-6b0c-4e55-8951-56f100c89842n%40googlegroups.com.

[ANTONIO LAUS]

Hi Michele,
First of all, thank you very much for your reply. I will immediately proceed to read the article in detail.
- From a first glance at the input files I seem to see some differences with the classical input files I used before, namely the presence of the FUNC entry (attached image). Is the entropic correction due to the funnel that usually ffs_tk does separately, or am I wrong?
- Also, I would like to ask whether in such a case it would be possible to investigate kinetics by varying the frequency (PACE), and if so, would the scripts used canonically for the other types of metadynamics be fine for calculating dwell times?

Best,

Antonio

cyl: 
DISTANCE
 
ATOMS
=v1,lig 
COMPONENTS
radius: 
MATHEVAL
 
ARG
=cyl.x,cyl.y 
FUNC
=sqrt(x*x+y*y 
PERIODIC
=NO

funnel: 
MATHEVAL
 
ARG
=radius,cyl.z 
VAR
=r,z 
FUNC
=(r+1.0*(-1.2+z))*step(-z+1.)+(r-0.2)*step(z-1 
PERIODIC
=NO
UPPER_WALLS
 
AT
=0 
ARG
=funnel 
KAPPA
=2000.0 
LABEL
=funnelwall  #funnel restraint
UPPER_WALLS
 
AT
=1.8 
ARG
=cyl.z 
KAPPA
=4000.0 
EXP
=2 
LABEL
=upper_wall  #upper limit of cyl.z

To view this discussion on the web visit https://groups.google.com/d/msgid/plumed-users/CAGBhpieR0i5_1eTqAbYDaR9mzTME34%2Bvjdx597DE0arroE3joA%40mail.gmail.com.

--

Antonio Laus, MSc

Fellow Researcher "Computer Aided Drug Design" 

Molecular Modelling & Drug Design Lab

Department of Life Sciences - University of Modena and Reggio Emilia

Via Giuseppe Campi 103

[Michele Invernizzi]

Hi Antonio,

Sorry, but I am not familiar with the `ffs_tk` tool you are mentioning. The MATHEVAL action (which is just a deprecated alias for CUSTOM) allows you to define a custom CV, defined by the given function FUNC. So that input is explicitly defining the funnel bias potential.

Regarding computing rates, it is in principle possible to use OPES in a similar way as infrequent metadynamics, but I suggest you have a look also at this approach https://pubs.acs.org/doi/abs/10.1021/acs.jctc.2c00806 It has been used also for ligand binging, for example https://www.nature.com/articles/s41467-022-33104-3
These papers share their input files, and some post-processing scripts that might be helpful. But you can also directly contact the authors, I'm sure they will be happy to help you out.

Best,

Michele

To view this discussion on the web visit https://groups.google.com/d/msgid/plumed-users/CAGZA7UnT2ZWd6nx%2BvCPAuRMvD72Evo9E_Hv5hWB5T8BQmdp%3D4w%40mail.gmail.com.

[ANTONIO LAUS]

Hi, I have read the articles and downloaded the files to understand more and for the calculation of residence time there should be no problem with the script, but I can't find the file for FES reweighting in the case of an OPES metadynamics in the presence of the funnel.
I would like to try to figure out if there is something I am missing:
- I assume that following the PLUMED masterclasses the reweighting of FES, for a standard WT-metadynamics with funnel, is done through a tool that they allow you to download (tooltip.tcl, funnel.tcl, ffs.tcl). So I do not understand if, in the absence of this tool, it can be done through a standard driver.dat or it would take a driver.dat with modifications including entropic correction of the funnel potential.
- Given this assumption, I am also doubtful in the case of OPES. That is, if when I go to do the driver counterpart via the reweighting script that is used for OPES, I would have to modify something to include the entropic correction of the funnel potential.

So unless it is somehow already calculated (during run), I would like to figure out how to do that.
I am attaching .dat files that I am using to do some testing, maybe you could help me .

#Funnel_Metadynamic

WHOLEMOLECULES ENTITY0=5,21,36,50,70,80,100,117,127,142,161,171,188,207,224,235,254,273,292,306,320,340,361,372,386,408,423,442,462,481,505,520,539,558,569,583,594,605,617,627,646,658,680,699,723,744,759,770,789,803,823,829,840,862,881,893,904,911,933,948,967,984,1003,1017,1036,1063,1069,1083,1105,1122,1134,1158,1172,1191,1205,1224,1240,1252,1266,1273,1292,1299,1316,1330,1352,1362,1374,1393,1412,1426,1440,1459,1466,1480,1499,1509,1531,1542,1549,1563,1585,1595,1615,1632,1647,1657,1676,1693,1703,1710,1720,1732,1751,1762,1779,1798,1805,1822,1842,1849,1865,1872,1892,1913,1924,1934,1955,1974,1990,2000,2015,2037,2053,2067,2083,2102,2116,2138,2155,2169,2181,2193,2208,2225,2246,2256,2280,2295,2306,2317,2327,2334,2341,2352,2372,2386,2402,2426,2440,2452,2466,2473,2496,2502,2519,2526,2550,2557,2571,2593,2609,2628,2647,2664,2683,2705,2720,2732,2749,2763,2778,2799,2818,2833,2848,2872,2896,2915,2937,2952,2971,2987,3009,3031,3048,3059,3076,3096,3115,3122,3151,3157,3176,3190,3209,3229,3245

t1: TORSION ATOMS=3259,3264,3265,3280
d1: DISTANCE ATOMS=1245,3263
lig: COM ATOMS=3259,3260,3261,3262,3263,3264,3265,3266,3267,3268,3269,3270,3271,3272,3273,3274,3275,3276,3277,3278,3279,3280,3281,3282,3283,3284,3285,3286,3287,3288,3289,3290,3291,3292,3293,3294,3295,3296,3297,3298,3299,3300,3301,3302,3303,3304,3305,3306,3307,3308,3309,3310

fps: FUNNEL_PS LIGAND=lig REFERENCE=4efu_CA_alignment.pdb ANCHOR=1245 POINTS=4.057,2.358,2.51,4.587,3.281,2.868

rmsd: RMSD REFERENCE=4efu_CA_alignment.pdb TYPE=OPTIMAL

FUNNEL ARG=fps.lp,fps.ld ZCC=2.0 ALPHA=0.8 RCYL=0.15 MINS=0.0 MAXS=3.7 KAPPA=80000 NBINS=1000 NBINZ=1000 FILE=BIAS LABEL=funnel

METAD ARG=d1,t1 SIGMA=0.01,0.05 HEIGHT=0.5 PACE=500 TEMP=300 BIASFACTOR=20 GRID_MIN=0.0,-pi GRID_MAX=8.0,pi GRID_WFILE=grid_w.dat CALC_RCT RCT_USTRIDE=1 GRID_WSTRIDE=25000000 LABEL=metad

LOWER_WALLS ARG=fps.lp AT=0.1 KAPPA=80000 EXP=2 OFFSET=0 LABEL=lwall
UPPER_WALLS ARG=rmsd AT=0.3 KAPPA=100000 EXP=2 OFFSET=0 LABEL=uwall-rmsd
UPPER_WALLS ARG=fps.lp AT=3.4 KAPPA=80000 EXP=2 OFFSET=0 LABEL=uwall
UPPER_WALLS ARG=d1 AT=3.4 KAPPA=80000 EXP=2 OFFSET=0 LABEL=distwall

PRINT STRIDE=500 ARG=* FILE=COLVAR

#exp OPES:
WHOLEMOLECULES ENTITY0=5,21,36,50,70,80,100,117,127,142,161,171,188,207,224,235,254,273,292,306,320,340,361,372,386,408,423,442,462,481,505,520,539,558,569,583,594,605,617,627,646,658,680,699,723,744,759,770,789,803,823,829,840,862,881,893,904,911,933,948,967,984,1003,1017,1036,1063,1069,1083,1105,1122,1134,1158,1172,1191,1205,1224,1240,1252,1266,1273,1292,1299,1316,1330,1352,1362,1374,1393,1412,1426,1440,1459,1466,1480,1499,1509,1531,1542,1549,1563,1585,1595,1615,1632,1647,1657,1676,1693,1703,1710,1720,1732,1751,1762,1779,1798,1805,1822,1842,1849,1865,1872,1892,1913,1924,1934,1955,1974,1990,2000,2015,2037,2053,2067,2083,2102,2116,2138,2155,2169,2181,2193,2208,2225,2246,2256,2280,2295,2306,2317,2327,2334,2341,2352,2372,2386,2402,2426,2440,2452,2466,2473,2496,2502,2519,2526,2550,2557,2571,2593,2609,2628,2647,2664,2683,2705,2720,2732,2749,2763,2778,2799,2818,2833,2848,2872,2896,2915,2937,2952,2971,2987,3009,3031,3048,3059,3076,3096,3115,3122,3151,3157,3176,3190,3209,3229,3245

t1: TORSION ATOMS=3259,3264,3265,3280
d1: DISTANCE ATOMS=1245,3263
lig: COM ATOMS=3259,3260,3261,3262,3263,3264,3265,3266,3267,3268,3269,3270,3271,3272,3273,3274,3275,3276,3277,3278,3279,3280,3281,3282,3283,3284,3285,3286,3287,3288,3289,3290,3291,3292,3293,3294,3295,3296,3297,3298,3299,3300,3301,3302,3303,3304,3305,3306,3307,3308,3309,3310

fps: FUNNEL_PS LIGAND=lig REFERENCE=test_CA_aligned.pdb ANCHOR=1245 POINTS=4.057,2.358,2.51,4.587,3.281,2.868

rmsd: RMSD REFERENCE=test_CA_aligned.pdb TYPE=OPTIMAL

FUNNEL ARG=fps.lp,fps.ld ZCC=2.0 ALPHA=0.8 RCYL=0.15 MINS=0.0 MAXS=3.7 KAPPA=80000 NBINS=1000 NBINZ=1000 FILE=BIAS LABEL=funnel

opes: OPES_METAD_EXPLORE ...
  ARG=d1,t1
  PACE=500
  BARRIER=80
  STATE_WFILE=STATE
  STATE_WSTRIDE=500*1 #0
  STORE_STATES
...

LOWER_WALLS ARG=fps.lp AT=0.1 KAPPA=80000 EXP=2 OFFSET=0 LABEL=lwall
UPPER_WALLS ARG=rmsd AT=0.3 KAPPA=100000 EXP=2 OFFSET=0 LABEL=uwall-rmsd
UPPER_WALLS ARG=fps.lp AT=3.4 KAPPA=80000 EXP=2 OFFSET=0 LABEL=uwall
UPPER_WALLS ARG=d1 AT=3.4 KAPPA=80000 EXP=2 OFFSET=0 LABEL=distwall

PRINT FMT=%g STRIDE=500 ARG=* FILE=COLVAR

Thank you very much,
Antonio 

[Michele Invernizzi]

Dear Antonio,

The entropic correction (I guess you are referring to eq. 1 of https://www.pnas.org/doi/epdf/10.1073/pnas.1303186110) is not part of the standard reweighting procedure, since it's only needed for the specific case of a funnel restraint potential. As such, it is not computed by the reweighting script, nor by standard plumed driver actions, such as REWEIGHT_BIAS. 

I see two options: 1. find a way to feed the reweighted fes to the ffs_tk tool, instead of its default sum_hills, 2. write your own few lines script to compute the entropic correction. 

Best,

Michele

PS
of course if someone reading this uses funnel metad, probably he could give you a working script or point you to a better solution

To view this discussion on the web visit https://groups.google.com/d/msgid/plumed-users/0a98bf0c-6298-4c82-b585-94860d303950n%40googlegroups.com.

[ANTONIO LAUS]

Dear Michele,

Thank you very much for the suggestion. I had actually tried using the FES_from_STATE script and then using the calculated FES in the ffs_tk tool for correction and it worked !
Instead as for reweighting scripts, I am using FESfromSTATE (although I have not fully understood, in terms of validity, the differences with FESformREWEIGHTING). I was wondering if it would be possible to set a stride as with sum_hills.
Thanks again for your patience!

Best,

Antonio

To view this discussion on the web visit https://groups.google.com/d/msgid/plumed-users/CAGBhpidsYCt03He68cyv2EZRKjLP9pi6%2B2468bf82P5%2BULWHLg%40mail.gmail.com.

[Michele Invernizzi]

Hi Antonio,

The FES_from_STATE gives you the FES estimate obtained as a function of the bias used during the simulation, while the FES_from_REWEIGTHING performs reweighting (with the assumption that the bias used is quasi-static). In my experience, they are extremely similar for OPES_METAD, and differ in OPES_METAD_EXPLORE, just like the sum_hills and reweighting of a METAD run (see fig. 4 of https://pubs.acs.org/doi/10.1021/acs.jctc.2c00152)

There is not a `--stride` option in the FES_from_STATE.py script, but there are two workarounds:

1.  use the keywork STORE_STATES in OPES_METAD, so that all the intermediate states are saved to the STATE file and then use the option `--all_stored`, this was done here
2. If you already run the simulation without saving the needed states, you can use the KERNELS file and the STATE_from_KERNELS script. You will have to write a small script that copies to a temporary file only the kernels that you want to include in the given FES.

Best,

Michele

To view this discussion on the web visit https://groups.google.com/d/msgid/plumed-users/CAGZA7Uk-zp4rMhJYgcaGr2kkyF1iWSYxJZQGgHcWX7DpqeCuWw%40mail.gmail.com.

[ANTONIO LAUS]

Dear Michele,
Thank you very much for your help and patience, you have clarified many doubts for me.
I only have two doubts from your last reply:
- If they are different FES_from_STATE from FES_from_REWEIGTHING in OPES_Metad_explore which one is advisable to consider?
- Also taking into consideration again the Funnel_reweighting that I am going to do with the ffs_tk tool, it occurs to me that it is convenient to use there FES_from_STATE and then do the reweighting via tool (in case I can skip the initial part and go to the point where the BIAS becomes quasi-static), please correct me if I am wrong.

Thanks again.

Best,

Antonio

To view this discussion on the web visit https://groups.google.com/d/msgid/plumed-users/CAGBhpifdoU0zGKvtUGjxQW01ave8G5hL9%2BrunsHK39h_mmjkJQ%40mail.gmail.com.

[Michele Invernizzi]

Dear Antonio,

They are different but statistically compatible, as shown in the paper. I prefer to use reweighting and perform multiple independent simulations to estimate the uncertainty, rather than relying on block average.

Best,

Michele

To view this discussion on the web visit https://groups.google.com/d/msgid/plumed-users/CAGZA7U%3Dd9H98%3DCQYDiu_%2BEWJT_Zsa8utVD%3DY0NLzT4i1j7SwcA%40mail.gmail.com.