Examine source of large error bars in 0.10.0 release

[jchodera]

The benchmark data for 0.10.0 has a few compounds with very large error bars.

We should take a look at these and see if there is something weird with the atom mappings, or other discernible issue.

[zhang-ivy]

I've gone through the release notes and tried to pick out potential reasons why the error has increased for some mutations:

• Bug generated by changes upstream in the openmm package -- default method for calculating solvent padding changed, which resulted in smaller boxes. Fixed by adding more padding to the simulation box, it is now 1.1 nm. Openmm issue Address change in padding behavior for addSolvent openmm/openmm#3502. Perses issue E openmm.OpenMMException: NonbondedForce: The cutoff distance cannot be greater than half the periodic box size. #949 (Add more padding to fix issue #949 #953)
• Use names (instead of indices) for fetching the force components of a system, avoiding issues with force reordering upstream in openmm. Issue validate_endstate_energies_point failing on openmm nightly #993 (grab forces by name instead of hard-coded index #976 and Specifying energy components by name, using dictionaries. #1007)
• Use of unique names for force components in HybridTopologyFactory (Rename second CustomBondForce in HybridTopologyFactory system and update docstrings #1022).
• Introduce unsampled_endstates boolean option in input YAML file, which enables/disables creation of unsampled endstates with long-range sterics correction. Issue Allow unsampled endstates to be controlled by yaml #1033 (unsampled endstates long range correction as an input yaml option #1037).

Nonte of these seem like they would be causing the increase in error bars though...

[ijpulidos]

I run the benchmark simulations without the endstates generation (unsampled_endstates: False on the template yaml) and I don't see the large error bars. So, I guess it is related to the unsampled endstates generation and likely introduced in #1023.

[zhang-ivy]

Ah great find! Could you point me to the code you're using to compute the free energies (and errors) for the above plots?

And also could you retrieve the .nc files for one transformation (both phases, with and without unsampled endstates)? We may need take a closer look at how the free energies (and errors) are computed.

[zhang-ivy]

LInk to the code used for computing free energies: https://github.com/choderalab/perses/blob/main/perses/analysis/load_simulations.py