salilab / pmi

Python Modeling Interface
https://integrativemodeling.org/nightly/doc/ref/namespaceIMP_1_1pmi.html
11 stars 11 forks source link

Different behavior between different runs when adding SAXS restraint #262

Closed ajtritt closed 7 months ago

ajtritt commented 7 months ago

I'm trying to load the EtfA and EtfB chain from this PDB structure, and then add a SAXS restraint.

Here is the code I am using: ```python import IMP import IMP.pmi as pmi import IMP.pmi.topology as topo import IMP.pmi.dof as dof import IMP.pmi.restraints.stereochemistry as sc from IMP.pmi.tools import shuffle_configuration from IMP.pmi.restraints import saxs fasta_path = "rcsb_pdb_7KOE.fasta" pdb_path = "7koe.pdb" saxs_path = '1_EtfAB_cryoEM/EtfAB_BlikeState_calculatedSAXS.dat' mdl = IMP.Model() sys = topo.System(mdl) cryoEM_st = sys.create_state() # load once so we can make a map to avoid having to use a huge key when we index below sequences = topo.Sequences(fasta_path) name_map = {k:v for k, v in zip(sequences.sequences.keys(), ['A', 'B', 'C', 'D'])} # and then reload so we can use simpler keys sequences = topo.Sequences(fasta_path, name_map=name_map) # Build hierarchy mol_etfB = cryoEM_st.create_molecule("EtfB", chain_id="B", sequence=sequences['B']) struct_B = mol_etfB.add_structure(pdb_path, chain_id='B', offset=-285) mol_etfA = cryoEM_st.create_molecule("EtfA", chain_id="A", sequence=sequences['A'][10:292]) struct_A = mol_etfA.add_structure(pdb_path, chain_id='A', offset=1) mol_etfB.add_representation(mol_etfB.get_atomic_residues(), resolutions=[1]) mol_etfB.add_representation(mol_etfB.get_non_atomic_residues(), resolutions=[1]) mol_etfA.add_representation(mol_etfA.get_atomic_residues(), resolutions=[1]) root_hier = sys.build() # Create rigid bodies etf_dof = dof.DegreesOfFreedom(mdl) rbN = etf_dof.create_rigid_body(mol_etfA[0:236] | mol_etfB[0:187]) fbA = etf_dof.create_flexible_beads(mol_etfA[236:251]) fbB = etf_dof.create_flexible_beads(mol_etfB[187:213]) rbC = etf_dof.create_rigid_body(mol_etfA[251:282] | mol_etfB[213:337]) # Add SAXS restrataint saxs_rest = saxs.SAXSRestraint([mol_etfA, mol_etfB], saxs_path) ```

When I run this, I get different errors each time. Below is an example of some of the errors. I suspect that I am doing something wrong but also wonder if the variable errors and the segfault are a symptom of a deeper issue.

Click here for the output ```bash (base) lm:IMP_Workshop$ cd /Users/ajtritt/projects/taskforce5/IMP/EtfABCX_IMP (base) lm:EtfABCX_IMP$ conda activate imp (imp) lm:EtfABCX_IMP$ python test.py > /dev/null Traceback (most recent call last): File "/Users/ajtritt/projects/taskforce5/IMP/EtfABCX_IMP/test.py", line 46, in saxs_rest = saxs.SAXSRestraint([mol_etfA, mol_etfB], saxs_path) ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^ File "/Users/ajtritt/anaconda3/envs/imp/lib/python3.11/site-packages/IMP/pmi/restraints/saxs.py", line 104, in __init__ self.restraint = IMP.saxs.Restraint(self.particles, self.profile, ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^ File "/Users/ajtritt/anaconda3/envs/imp/lib/python3.11/site-packages/IMP/saxs/__init__.py", line 1260, in __init__ _IMP_saxs.Restraint_swiginit(self, _IMP_saxs.new_Restraint(*args)) ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^ _IMP_kernel.ValueException: Atom is not the child of a residue "N" (Un) (imp) lm:EtfABCX_IMP$ python test.py > /dev/null Traceback (most recent call last): File "/Users/ajtritt/projects/taskforce5/IMP/EtfABCX_IMP/test.py", line 46, in saxs_rest = saxs.SAXSRestraint([mol_etfA, mol_etfB], saxs_path) ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^ File "/Users/ajtritt/anaconda3/envs/imp/lib/python3.11/site-packages/IMP/pmi/restraints/saxs.py", line 104, in __init__ self.restraint = IMP.saxs.Restraint(self.particles, self.profile, ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^ File "/Users/ajtritt/anaconda3/envs/imp/lib/python3.11/site-packages/IMP/saxs/__init__.py", line 1260, in __init__ _IMP_saxs.Restraint_swiginit(self, _IMP_saxs.new_Restraint(*args)) ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^ _IMP_kernel.InternalException: Corrupted Key Table asking for key 1735289198 with a table of size 536 (imp) lm:EtfABCX_IMP$ python test.py > /dev/null Traceback (most recent call last): File "/Users/ajtritt/projects/taskforce5/IMP/EtfABCX_IMP/test.py", line 46, in saxs_rest = saxs.SAXSRestraint([mol_etfA, mol_etfB], saxs_path) ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^ File "/Users/ajtritt/anaconda3/envs/imp/lib/python3.11/site-packages/IMP/pmi/restraints/saxs.py", line 104, in __init__ self.restraint = IMP.saxs.Restraint(self.particles, self.profile, ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^ File "/Users/ajtritt/anaconda3/envs/imp/lib/python3.11/site-packages/IMP/saxs/__init__.py", line 1260, in __init__ _IMP_saxs.Restraint_swiginit(self, _IMP_saxs.new_Restraint(*args)) ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^ _IMP_kernel.InternalException: Corrupted Key Table asking for key 9234 with a table of size 536 (imp) lm:EtfABCX_IMP$ python test.py > /dev/null Segmentation fault: 11 ```

Links to inputs if you need them:

benmwebb commented 7 months ago

Short answer: change the last line to read

saxs_rest = saxs.SAXSRestraint([mol_etfA, mol_etfB], saxs_path, ff_type=IMP.saxs.RESIDUES)

Longer answer: the SAXS restraint by default works with atomic form factors, and you have provided (as is typical for PMI) only 1-bead-per-residue representation. Sure, you're doing something wrong, but our documentation and error messages should definitely do a better job here. I can't reproduce the table corruption or segfaults you're seeing, but almost certainly that's caused by the restraint dumbly taking the Residue particles you've provided and trying to treat them as Atom particles, without checking. I will address that.

ajtritt commented 7 months ago

@benmwebb Thanks for digging into this. I actually wondered if ff_type had something to do with it, but the options in the docs didn't make any sense. After setting the right form factor, everything seems to be running smoothly.