OpenFreeEnergy / openfe_analysis

GNU General Public License v3.0
2 stars 1 forks source link

Easy access to non-hybrid topologies for endstates #34

Open tlhr opened 1 week ago

tlhr commented 1 week ago

Is there an easy way to get access to a clean representation (non-hybrid) of the ligand trajectory at lambda[0] and lambda[-1]? Maybe I'm overcomplicating the indexing procedure, but actually it seems quite tedious.

In an ideal world I would envision something like an extra flag for the reader, e.g.:

u = mda.Universe('hybrid_system.pdb', 'simulation.nc', format=FEReader, state_id=0, hybrid=False)

This would make subsequent analyses a lot easier!

tlhr commented 1 week ago

I just realized that I was indeed overcomplicating the indexing, it seems that the lambda state is stored in the bfactor field (with 0.25, 0.5, 0.75 corresponding to state A, both states, and state B respectively). So getting a clean trajectory just involves something like this:

atom_mapping_b = np.array([atom.ix for atom in u.atoms if atom.bfactor in (0.5, 0.75)])
state_b = u.atoms[atom_mapping_b]

Leaving this open to see what the feelings on a flag are, otherwise this can be closed!

IAlibay commented 1 week ago

Ah yes, thanks for raising this @tlhr - the bfactor hack is.. well kind of a hack, we probably need to find a better long term solution to this problem and/or make sure it's well documented!