Early alchemiscale results

[IAlibay]

Since @dotsdl has both a poster at the free energy workshop and a talk at the OMSF workshop, we thought it would be good to aggregate initial results here.

[IAlibay]

Update on my end: we should have a full set of openff 2.0 and 2.1-rc1 comparisons for TYK2, Thrombin, and MCL-1 by ~ Thursday (I'm flying tomorrow so won't be able to aggregate).

Input data:

• For an initial "backwards compatibility" check, we're using the data from https://github.com/OpenFreeEnergy/openfe-benchmarks/tree/main/openfe_benchmarks/data (future iterations will use the PLB data directly).

Protocol details:

• 3 replicas
• HREX w/ 11 windows
• 2 ns equilibration and "up to" 5 ns production (early termination once the error drops below 0.2 kT)
• HMR (3 amu) + 4 fs timestep
• cutoff is set at 0.9 nm, we aren't doing unsampled edges

[IAlibay]

All results

[dotsdl]

Awesome! Thank you for this @IAlibay!!!

[dotsdl]

@IAlibay any thoughts on the outliers on p38? in the DDG plot (right)?

[IAlibay]

@dotsdl see here: https://openforcefieldgroup.slack.com/archives/CKSHCE7SB/p1683813489278679?thread_ts=1683813489.278679&cid=CKSHCE7SB

Short answer is that we don't know right now (I don't have tools to know what the FF changes were for these two ligands).

[dotsdl]

Ah cool, just wanted to make sure I was able to answer the question if I get it tomorrow; much appreciated!

[ijpulidos]

These are the first results I obtained running the NEQ cycling protocol from perses. I'm uploading the figures and the arsenic/cinnabar csv file that I used for them.

These are results from the tyk2 system in https://github.com/OpenFreeEnergy/openfe-benchmarks/tree/main/openfe_benchmarks/data

tyk2_arsenic.csv

[IAlibay]

Ok, here's a comparison between HREX and NECycling (re-plotted because I changed the plotting style and I didn't want that to be the source of confusion):

The HREX results here use the same protocol details as those outlined here: https://github.com/orgs/openforcefield/discussions/24#discussioncomment-5722528

Note the 0.2 kT label represents the early termination threshold (i.e. it'll run for up to 11 ns per window, or until a 0.2 kT threshold is reached, whatever happens first).

It looks like all the problems seem to stem from just one edge: lig_ejm_31 -> lig_ejm_54, which is overestimated by the NECycling protocol.

Just to check that this wasn't a sampling issue, I reran the TYK2 system without turning on early termination (labelled as "11 * 7 ns"). Those are consistent.

We could try running the repex protocol much longer for the one edge, maybe it's some kind of slow loop motion?

[IAlibay]

@ijpulidos - is it possible that this used the old gufe default for hmass? I've seen bad edges on 4 amu mass HMR just randomly pop up with relative hydration free energies.

[dotsdl]

We could try running the repex protocol much longer for the one edge, maybe it's some kind of slow loop motion?

@IAlibay did you mean here that we should try running longer with the neqcyc protocol, not repex, since we see the overestimation in the neqcyc result for that one edge?

[IAlibay]

@dotsdl - yeah sorry neqcyc should be redone in the first instance. But we should probably run both to exhaustion just to check tbh.

[ijpulidos]

@ijpulidos - is it possible that this used the old gufe default for hmass? I've seen bad edges on 4 amu mass HMR just randomly pop up with relative hydration free energies.

I did use gufe 0.7.2 which is not the latest release. But I'm explicitely setting settings.forcefield_settings.hydrogen_mass = 3.0, so it should be using the correct value.

@dotsdl - yeah sorry neqcyc should be redone in the first instance. But we should probably run both to exhaustion just to check tbh.

For now I think a quick way to check is re-running the "problematic" lig_ejm_31 -> lig_ejm_54 transformation for a longer time, I will double check the performance we are getting and increase the simulation times accordingly.

[ppxasjsm]

Wow these are really cool!