Hello,
I am modeling a multistate system in PMI using CrossLinkingMassSpectrometryRestraint in imp-2.15, and trying to adapt imp-sampcon for the purpose.
I set up my system with 2 states having identical components with:
topology = IMP.pmi.topology.TopologyReader(topology_file, fasta_dir=datadirectory, pdb_dir=datadirectory)
topology2 = IMP.pmi.topology.TopologyReader(topology_file, fasta_dir=datadirectory, pdb_dir=datadirectory)
bm = IMP.pmi.macros.BuildSystem(m) bm.add_state(topology) bm.add_state(topology2)
root_hierarchy, dof = bm.execute_macro(max_rb_trans=rb_max_trans, max_rb_rot=rb_max_rot, max_bead_trans=bead_max_trans)
and later set up crosslinking restraint (after creating crosslink database) by:
cl_dsso = IMP.pmi.restraints.crosslinking.CrossLinkingMassSpectrometryRestraint(root_hier=root_hierarchy,
database=cldb2,
length=21.0,
slope=0.01,
resolution=1.0,
filelabel='dsso',
label="dsso",
linker=ihm.cross_linkers.dsso,
weight=1)
After running replica exchange, i am trying to extract distances for my crosslinks (in theory 2 per protein copy, one in state 0 and one in state 1). However, when I for example extract distances from the stat file by
stat_file=IMP.pmi.output.ProcessOutput('rmfs/0.rmf3') #extract distances for a specific crosslink restraint stat_file.get_fields(['CrossLinkingMassSpectrometryRestraint_Distance_|dsso|1278.1|Protein1|148|Protein2|689|1|PSI|'])
I get a single distance per frame, corresponding to the distance in the latest state to be added, not to the shortest distance in the ensemble, and not one distance per state, as would make sense for analysing crosslink violations. As far as I can tell the score is behaving correctly overall (optimising for one state to be below the cutoff).
I can still successfully run the first part of sampcon using score cutoffs rather than distances. But I wonder if it is possible to extract the shortest distance for each crosslink? Have I set up the multistate system incorrectly? (I am trying to work out how to set up the alignment etc for rmsd in such a case, but I have a long way to go there.)
I can provide a test case if needed.
Many thanks in advance,
Andrea