Passing a hierarchy to stereochemistry.ExcludedVolumeSphere does not work! (IMP 2.7.0)

[fergaljd]

I've been running my simulations using pmi, and the IMP.pmi.macros.BuildSystem class before adding restraints, and I've noticed that passing a hierarchy to the ExcludedVolumeSphere object will fail silently.

model = IMP.Model()
bs = IMP.pmi.macros.BuildSystem(model)
bs.add_state("topology.txt")
root_hier, tfiid_dof = bs.execute_macro()
#Below doesn't work!
ev = IMP.pmi.stereochemistry.ExcludedVolumeSphere(included_objects=root_hier)
ev.add_to_model()
ev.evaluate() # Will always be 0.0!

I think this behaviour comes from the following calls in the ExcludedVolumeSphere constructor

hierarchies = IMP.pmi.tools.input_adaptor(included_objects, resolution, flatten=True)
...
included_ps = [h.get_particle() for h in hierarchies]

Passing the IMP.pmi.tools.input_adaptor function an IMP.atom.hierarchy hierarchy just returns the hierarchy itself in a list, which then leads to only a single particle getting set in the restraint. I've worked around it for now by using:

ev = stereochemistry.ExcludedVolumeSphere(included_objects = bs.system)

but I think the class should either work properly with hierarchies, or not accept them at all

[Pellarin]

That should have been posted in the IMP.pmi issues. I'll open an issue for you . Can you please define what you mean by "it is not working?". Anyways, your question is interesting, because apparently IMP.pmi.tools.input_adaptor does not return the leaves of the input hierarchies.

Related to https://github.com/salilab/pmi/issues/234

[fergaljd]

I mean that if you pass the root hierarchy to IMP.pmi.stereochemistry.ExcludedVolumeSphere, the restraint will always evaluate to 0 during simulation. Since IMP.pmi.tools.input_adaptor simply returns the input hierarchy, only a single particle ID, representing the System is used to search for volume overlap.

I'm guessing it is due to the following section the input_adapter() function in IMP.pmi.tools.py

 ...
 hier_list = []
 if tp in (IMP.pmi.topology.System,IMP.pmi.topology.State,
           IMP.pmi.topology.Molecule,IMP.pmi.topology.TempResidue):
     # if PMI, perform selection using gathered indexes
     <snip>

 else:
     try:
         if IMP.atom.Hierarchy.get_is_setup(stuff[0]):
             if flatten:
                 hier_list = stuff
             else:
                 hier_list = [stuff]

[Pellarin]

OK, you get 0 score, where you expect to have non-zero.

[cgreenberg]

Right, well, this is a feature not a bug :)

The idea here is that passing PMI objects is the way to use PMI restraints and other classes. So if you pass a PMI System, State, or Molecule, or a collection of them, the restraint will in this case get all components particles and create an excluded volume restraint on them. This is good for short setup scripts.

However, in some cases advanced users may want to customize exactly which particles are assigned the restraint. In that case you pass IMP hierarchies or particles. e.g. you could use IMP.atom.Selection to get a specific chain, residue type, resolution, etc etc. In this case the restraint will be "dumb" and only create the restraint on the exact passed hierarchies. So it's up to you to select the particles first.

Clearly we need to document this better but the dual-mode easy/hard way of using PMI restraints is intentional.

[cgreenberg]

We should probably change the return value for BuildSystem.execute_macro() to give the PMI System instead of the root hierarchy though.

[Pellarin]

OK, for the present case:

hiers=IMP.atom.Selection(root_hier,resolution=10).get_selected_particles()
ev = IMP.pmi.stereochemistry.ExcludedVolumeSphere(included_objects=hiers)

should have returned the same result as passing bs.system.

1. Having input_adaptor realeasing the leaves of input hierarchies should not conflict with present code.

2. I'm starting to be more keen to use only IMP.atom.Hierarchies, and less IMP.pmi.topology objects.