New matching function for convenience superposition: MatchResidueByIdx

modules/mol/alg/doc/molalg.rst

@@ -46,6 +46,8 @@
 
 .. autofunction:: MatchResidueByNum
 
+.. autofunction:: MatchResidueByIdx
+
 .. autofunction:: Superpose
 
 

modules/mol/alg/pymod/superpose.py

@@ -27,6 +27,7 @@ def ParseAtomNames(atoms):
   :param atoms: Identifier or list of atoms
   :type atoms: :class:`str`, :class:`list`, :class:`set`
   """
+  ## get a set of atoms or None
   if atoms==None:
     return None
   if isinstance(atoms, str):
@@ -62,6 +63,16 @@ def _fetch_atoms(r_a, r_b, result_a, result_b, atmset):
   return result_a, result_b
 
 
+def _no_of_chains(ent_a, ent_b):
+  """
+  for internal use, only
+  """
+  ## get lower no. of chains
+  if ent_a.chain_count < ent_b.chain_count:
+    return ent_a.chain_count
+  return ent_b.chain_count
+
+
 def MatchResidueByNum(ent_a, ent_b, atoms='all'):
   """
   Returns a tuple of views containing exactly the same number of atoms.
@@ -82,20 +93,15 @@ def MatchResidueByNum(ent_a, ent_b, atoms='all'):
   ## init. final views
   result_a=_EmptyView(ent_a)
   result_b=_EmptyView(ent_b)
-  ## get lower no. of chains
-  if ent_a.chain_count < ent_b.chain_count:
-    n_chains=ent_a.chain_count
-  else:
-    n_chains=ent_b.chain_count
-  ## get a set of atoms or None
+  n_chains=_no_of_chains(ent_a, ent_b)
   atmset=ParseAtomNames(atoms)
   ## iterate chains
   for i in range(0, n_chains):
     chain_a=ent_a.chains[i]
     chain_b=ent_b.chains[i]
-    ## residues of chains need to be consecutively numbered (sorted)
+    residues_a=iter(chain_a.residues)
+    ## decide on order of residues
     if chain_a.InSequence() and chain_b.InSequence():
-      residues_a=iter(chain_a.residues)
       residues_b=iter(chain_b.residues)
       ## check residues & copy to views
       try:
@@ -112,7 +118,6 @@ def MatchResidueByNum(ent_a, ent_b, atoms='all'):
         pass
     else:
       ## iterate one list of residues, search in other list
-      residues_a=iter(chain_a.residues)
       try:
         while True:
           r_a=residues_a.next()
@@ -126,6 +131,46 @@ def MatchResidueByNum(ent_a, ent_b, atoms='all'):
   return result_a, result_b
 
 
+def MatchResidueByIdx(ent_a, ent_b, atoms='all'):
+  """
+  Returns a tuple of views containing exactly the same number of atoms.
+  Residues are matched by position in the chains of an entity. A subset of
+  atoms to be included in the views can be specified in the **atoms** argument.
+  Regardless of what the list of **atoms** says, only those present in two
+  matched residues will be included in the views. Chains are processed in order
+  of appearance. If **ent_a** and **ent_b** contain a different number of
+  chains, processing stops with the lower count. The number of residues per
+  chain is supposed to be the same.
+
+  :param ent_a: The first entity
+  :type ent_a: :class:`~ost.mol.EntityView` or :class:`~ost.mol.EntityHandle`
+  :param ent_b: The second entity
+  :type ent_b: :class:`~ost.mol.EntityView` or :class:`~ost.mol.EntityHandle`
+  :param atoms: The subset of atoms to be included in the two views.
+  :type atoms: :class:`str`, :class:`list`, :class:`set`
+  """
+  not_supported="MatchResidueByIdx has no support for chains of different "\
+               +"lengths"
+  ## init. final views
+  result_a=_EmptyView(ent_a)
+  result_b=_EmptyView(ent_b)
+  n_chains=_no_of_chains(ent_a, ent_b)
+  atmset=ParseAtomNames(atoms)
+  ## iterate chains
+  for i in range(0, n_chains):
+    chain_a=ent_a.chains[i]
+    chain_b=ent_b.chains[i]
+    ## check equal no. of residues
+    if chain_a.residue_count!=chain_b.residue_count:
+      raise RuntimeError(not_supported)
+    ## iterate residues & copy to views
+    for r_a,r_b in zip(chain_a.residues, chain_b.residues):
+      result_a,result_b=_fetch_atoms(r_a, r_b, result_a, result_b, atmset)
+  result_a.AddAllInclusiveBonds()
+  result_b.AddAllInclusiveBonds()
+  return result_a, result_b
+
+
 def Superpose(ent_a, ent_b, match='number', atoms='all'):
   """
   Superposes the first entity onto the second. To do so, two views are created,
@@ -137,6 +182,9 @@ def Superpose(ent_a, ent_b, match='number', atoms='all'):
   * ``number`` - select residues by residue number, includes **atoms**, calls
     :func:`~ost.mol.alg.MatchResidueByNum`
 
+  * ``index`` - select residues by index in chain, includes **atoms**, calls
+    :func:`~ost.mol.alg.MatchResidueByIdx`
+
   :param ent_a: The first entity
   :type ent_a: :class:`~ost.mol.EntityView` or :class:`~ost.mol.EntityHandle`
   :param ent_b: The second entity
@@ -150,6 +198,8 @@ def Superpose(ent_a, ent_b, match='number', atoms='all'):
   ## create views to superpose
   if match.upper()=='NUMBER':
     view_a, view_b=MatchResidueByNum(ent_a, ent_b, atoms)
+  elif match.upper()=='INDEX':
+    view_a, view_b=MatchResidueByIdx(ent_a, ent_b, atoms)
   else:
     raise ValueError(not_supported)
   ## action

modules/mol/alg/tests/test_convenient_superpose.py

@@ -7,11 +7,6 @@ class TestConvenientSuperpose(unittest.TestCase):
   def setUp(self):
     pass
   
-  def runAtomOrdering(self, view1, view2):
-    # call atom ordering function here
-    view1, view2 = ost.mol.alg.MatchResidueByNum(view1, view2)
-    return view1, view2
-  
   def assertEqualAtomOrder(self, view1, view2):
     self.assertEquals(len(view1.atoms),len(view2.atoms))
     for a1, a2 in zip(view1.atoms, view2.atoms):
@@ -30,45 +25,74 @@ class TestConvenientSuperpose(unittest.TestCase):
   def testCorrectlyOrdered(self):
     ent1_ent = io.LoadEntity(os.path.join("testfiles","1aho.pdb"))
     ent1_view = ent1_ent.Select("")
-    view1, view2 = self.runAtomOrdering(ent1_ent, ent1_ent)
+    ## test MatchResidueByNum
+    view1, view2 = mol.alg.MatchResidueByNum(ent1_ent, ent1_ent)
+    self.assertEqualAtomOrder(view1, view2)    
+    view1, view2 = mol.alg.MatchResidueByNum(ent1_view, ent1_ent)
+    self.assertEqualAtomOrder(view1, view2)
+    view1, view2 = mol.alg.MatchResidueByNum(ent1_ent, ent1_view)
+    self.assertEqualAtomOrder(view1, view2)
+    view1, view2 = mol.alg.MatchResidueByNum(ent1_view, ent1_view)
+    self.assertEqualAtomOrder(view1, view2)
+    ## test MatchResidueByIdx
+    view1, view2 = mol.alg.MatchResidueByIdx(ent1_ent, ent1_ent)
     self.assertEqualAtomOrder(view1, view2)    
-    view1, view2 = self.runAtomOrdering(ent1_view, ent1_ent)
+    view1, view2 = mol.alg.MatchResidueByIdx(ent1_view, ent1_ent)
     self.assertEqualAtomOrder(view1, view2)
-    view1, view2 = self.runAtomOrdering(ent1_ent, ent1_view)
+    view1, view2 = mol.alg.MatchResidueByIdx(ent1_ent, ent1_view)
     self.assertEqualAtomOrder(view1, view2)
-    view1, view2 = self.runAtomOrdering(ent1_view, ent1_view)
+    view1, view2 = mol.alg.MatchResidueByIdx(ent1_view, ent1_view)
     self.assertEqualAtomOrder(view1, view2)
     
   def testMissingFirstAtom(self):
     ent_view = io.LoadEntity(os.path.join("testfiles","1aho.pdb")).Select("")
     ent_view_missing = ent_view.Select("not (cname=A and rnum=1 and aname=N)")
-    view1, view2 = self.runAtomOrdering(ent_view, ent_view_missing)
+    ## test MatchResidueByNum
+    view1, view2 = mol.alg.MatchResidueByNum(ent_view, ent_view_missing)
+    self.assertEqualAtomOrder(view1, view2)
+    view1, view2 = mol.alg.MatchResidueByNum(ent_view_missing, ent_view)
     self.assertEqualAtomOrder(view1, view2)
-    view1, view2 = self.runAtomOrdering(ent_view_missing, ent_view)
+    ## test MatchResidueByIdx
+    view1, view2 = mol.alg.MatchResidueByIdx(ent_view, ent_view_missing)
+    self.assertEqualAtomOrder(view1, view2)
+    view1, view2 = mol.alg.MatchResidueByIdx(ent_view_missing, ent_view)
     self.assertEqualAtomOrder(view1, view2)
     
   def testMissingManyAtoms(self):
     ent_view = io.LoadEntity(os.path.join("testfiles","1aho.pdb")).Select("")
     ent_view_missing = ent_view.Select("not (cname=A and rnum=3,19,32 and aname=CB,CA,CD)")
-    view1, view2 = self.runAtomOrdering(ent_view, ent_view_missing)
+    ## test MatchResidueByNum
+    view1, view2 = mol.alg.MatchResidueByNum(ent_view, ent_view_missing)
+    self.assertEqualAtomOrder(view1, view2)
+    view1, view2 = mol.alg.MatchResidueByNum(ent_view_missing, ent_view)
     self.assertEqualAtomOrder(view1, view2)
-    view1, view2 = self.runAtomOrdering(ent_view_missing, ent_view)
+    ## test MatchResidueByIdx
+    view1, view2 = mol.alg.MatchResidueByIdx(ent_view, ent_view_missing)
+    self.assertEqualAtomOrder(view1, view2)
+    view1, view2 = mol.alg.MatchResidueByIdx(ent_view_missing, ent_view)
     self.assertEqualAtomOrder(view1, view2)
     
   def testMissingFirstResidue(self):
     ent_view = io.LoadEntity(os.path.join("testfiles","1aho.pdb")).Select("")
     ent_view_missing = ent_view.Select("not (cname=A and rnum=1)")
-    view1, view2 = self.runAtomOrdering(ent_view, ent_view_missing)
+    ## test MatchResidueByNum
+    view1, view2 = mol.alg.MatchResidueByNum(ent_view, ent_view_missing)
     self.assertEqualAtomOrder(view1, view2)
-    view1, view2 = self.runAtomOrdering(ent_view_missing, ent_view)
+    view1, view2 = mol.alg.MatchResidueByNum(ent_view_missing, ent_view)
     self.assertEqualAtomOrder(view1, view2)
-  
+
   def testMissingHydrogens(self):
     ent_view = io.LoadEntity(os.path.join("testfiles","1aho.pdb")).Select("")
     ent_view_missing = ent_view.Select("ele!=H")
-    view1, view2 = self.runAtomOrdering(ent_view, ent_view_missing)
+    ## test MatchResidueByNum
+    view1, view2 = mol.alg.MatchResidueByNum(ent_view, ent_view_missing)
+    self.assertEqualAtomOrder(view1, view2)
+    view1, view2 = mol.alg.MatchResidueByNum(ent_view_missing, ent_view)
     self.assertEqualAtomOrder(view1, view2)
-    view1, view2 = self.runAtomOrdering(ent_view_missing, ent_view)
+    ## test MatchResidueByIdx
+    view1, view2 = mol.alg.MatchResidueByIdx(ent_view, ent_view_missing)
+    self.assertEqualAtomOrder(view1, view2)
+    view1, view2 = mol.alg.MatchResidueByIdx(ent_view_missing, ent_view)
     self.assertEqualAtomOrder(view1, view2)
   
   def testWrongAtomOrder(self):
@@ -82,9 +106,15 @@ class TestConvenientSuperpose(unittest.TestCase):
         random.shuffle(atoms)
         for a in atoms:
           ent_view_wrong.AddAtom(a)
-    view1, view2 = self.runAtomOrdering(ent_view, ent_view_wrong)
+    ## test MatchResidueByNum
+    view1, view2 = mol.alg.MatchResidueByNum(ent_view, ent_view_wrong)
+    self.assertEqualAtomOrder(view1, view2)
+    view1, view2 = mol.alg.MatchResidueByNum(ent_view_wrong, ent_view)
+    self.assertEqualAtomOrder(view1, view2)
+    ## test MatchResidueByIdx
+    view1, view2 = mol.alg.MatchResidueByIdx(ent_view, ent_view_wrong)
     self.assertEqualAtomOrder(view1, view2)
-    view1, view2 = self.runAtomOrdering(ent_view_wrong, ent_view)
+    view1, view2 = mol.alg.MatchResidueByIdx(ent_view_wrong, ent_view)
     self.assertEqualAtomOrder(view1, view2)
   
   def testWrongResidueOrder(self): 
@@ -98,9 +128,10 @@ class TestConvenientSuperpose(unittest.TestCase):
         ent_view_wrong.AddResidue(r)
         for a in r.atoms:
           av=ent_view_wrong.AddAtom(a)
-    view1, view2 = self.runAtomOrdering(ent_view, ent_view_wrong)
+    ## test MatchResidueByNum
+    view1, view2 = mol.alg.MatchResidueByNum(ent_view, ent_view_wrong)
     self.assertEqualAtomOrder(view1, view2)
-    view1, view2 = self.runAtomOrdering(ent_view_wrong, ent_view)
+    view1, view2 = mol.alg.MatchResidueByNum(ent_view_wrong, ent_view)
     self.assertEqualAtomOrder(view1, view2)
     
 if __name__ == "__main__":