From 0797fe2956168479fb9ea90149bbf5cb15687f66 Mon Sep 17 00:00:00 2001
From: Gerardo Tauriello <gerardo.tauriello@unibas.ch>
Date: Mon, 9 May 2016 18:19:48 +0200
Subject: [PATCH] SCHWED-877: add special treatment for peptide ligands
---
sidechain/pymod/_reconstruct_sidechains.py | 126 ++++++++++++++-------
1 file changed, 82 insertions(+), 44 deletions(-)
diff --git a/sidechain/pymod/_reconstruct_sidechains.py b/sidechain/pymod/_reconstruct_sidechains.py
index bb8f5410..7f8cbabe 100644
--- a/sidechain/pymod/_reconstruct_sidechains.py
+++ b/sidechain/pymod/_reconstruct_sidechains.py
@@ -15,6 +15,46 @@ def _GetRotamerIDs(res_list):
rotamer_ids[i] = rot_id
return rotamer_ids
+def _GetPhiAngle(r):
+ '''Extract phi angle for residue r.'''
+ # def. fallback = helix
+ phi = -1.0472
+ # try to get phi from torsion angles
+ tor = r.GetPhiTorsion()
+ if tor.IsValid():
+ phi = tor.GetAngle()
+ else:
+ r_prev = r.handle.prev
+ if r_prev.IsValid() and mol.InSequence(r_prev, r.handle):
+ c_prev = r_prev.FindAtom("C")
+ n = r.FindAtom("N")
+ ca = r.FindAtom("CA")
+ c = r.FindAtom("C")
+ if c_prev.IsValid() and n.IsValid() and ca.IsValid() and c.IsValid():
+ phi = geom.DihedralAngle(c_prev.GetPos(),n.GetPos(),
+ ca.GetPos(),c.GetPos())
+ return phi
+
+def _GetPsiAngle(r):
+ '''Extract psi angle for residue r.'''
+ # def. fallback = helix
+ psi = -0.7854
+ # try to get psi from torsion angles
+ tor = r.GetPsiTorsion()
+ if tor.IsValid():
+ psi = tor.GetAngle()
+ else:
+ r_next = r.handle.next
+ if r_next.IsValid() and mol.InSequence(r.handle, r_next):
+ n = r.FindAtom("N")
+ ca = r.FindAtom("CA")
+ c = r.FindAtom("C")
+ n_next = r_next.FindAtom("N")
+ if n.IsValid() and ca.IsValid() and c.IsValid() and n_next.IsValid():
+ psi = geom.DihedralAngle(n.GetPos(), ca.GetPos(),
+ c.GetPos(), n_next.GetPos())
+ return psi
+
def _GetDihedrals(res_list):
'''Extract dihedral angles for all residues.
Returns phi and psi angles as 2 lists with length = len(res_list).
@@ -22,42 +62,8 @@ def _GetDihedrals(res_list):
phi_angles = [0.0] * len(res_list)
psi_angles = [0.0] * len(res_list)
for i,r in enumerate(res_list):
- # def. fallback = helix
- phi = -1.0472
- psi = -0.7854
-
- # try to get phi from torsion angles
- tor = r.GetPhiTorsion()
- if tor.IsValid():
- phi = tor.GetAngle()
- elif i > 0:
- r_prev = res_list[i-1]
- if mol.InSequence(r_prev.handle,r.handle):
- c_prev = r_prev.FindAtom("C")
- n = r.FindAtom("N")
- ca = r.FindAtom("CA")
- c = r.FindAtom("C")
- if c_prev.IsValid() and n.IsValid() and ca.IsValid() and c.IsValid():
- phi = geom.DihedralAngle(c_prev.GetPos(),n.GetPos(),
- ca.GetPos(),c.GetPos())
-
- # try to get psi from torsion angles
- tor = r.GetPsiTorsion()
- if tor.IsValid():
- psi = tor.GetAngle()
- elif i < (len(res_list) - 1):
- r_next = res_list[i+1]
- if mol.InSequence(r.handle,r_next.handle):
- n = r.FindAtom("N")
- ca = r.FindAtom("CA")
- c = r.FindAtom("C")
- n_next = r_next.FindAtom("N")
- if n.IsValid() and ca.IsValid() and c.IsValid() and n_next.IsValid():
- psi = geom.DihedralAngle(n.GetPos(), ca.GetPos(),
- c.GetPos(), n_next.GetPos())
- # store
- phi_angles[i] = phi
- psi_angles[i] = psi
+ phi_angles[i] = _GetPhiAngle(r)
+ psi_angles[i] = _GetPsiAngle(r)
return phi_angles, psi_angles
def _AddBackboneFrameResidues(frame_residues, res_list, rotamer_ids,
@@ -76,14 +82,46 @@ def _AddLigandFrameResidues(frame_residues, ent_lig, rotamer_settings, offset):
'''Update frame_residues (list) with FrameResidues for res. in ent_lig.
Set offset >= number of non-ligand residues (used for residue_index).
'''
- # TODO? Split peptides and non-peptides?
- # for now: add as static carbons
- for i,r in enumerate(ent_lig.residues):
- try:
- frame_residue = sidechain.ConstructFrameResidue(r.handle, offset+i)
- frame_residues.append(frame_residue)
- except:
- continue
+ # parse ligand residues
+ for i, res in enumerate(ent_lig.residues):
+ res_idx = offset + i
+ is_done = False
+ # special treatment for peptides
+ if res.IsPeptideLinking():
+ rot_id = sidechain.TLCToRotID(res.GetName())
+ if rot_id != sidechain.XXX:
+ # get more info
+ phi = _GetPhiAngle(res)
+ r_prev = res.handle.prev
+ n_ter = not r_prev.IsValid() \
+ or not mol.InSequence(r_prev, res.handle)
+ r_next = res.handle.next
+ c_ter = not r_next.IsValid() \
+ or not mol.InSequence(res.handle, r_next)
+ # try to add frame residues (ignore exceptions)
+ try:
+ fr1 = sidechain.ConstructBackboneFrameResidue(\
+ res.handle, rot_id, res_idx, rotamer_settings,
+ phi, n_ter, c_ter)
+ if rot_id != sidechain.ALA and rot_id != sidechain.GLY:
+ fr2 = sidechain.ConstructSidechainFrameResidue(\
+ res.handle, rot_id, res_idx, rotamer_settings)
+ frame_residues.extend([fr1,fr2])
+ else:
+ frame_residues.append(fr1)
+ except:
+ pass # ignore peptide treatment and treat below
+ else:
+ is_done = True
+ # if it failed, treat it as an unknown entity
+ if not is_done:
+ # try to add frame residues (skip exceptions)
+ try:
+ # TODO: allow non-carbon treatment
+ fr = sidechain.ConstructFrameResidue(res.handle, res_idx)
+ frame_residues.append(fr)
+ except:
+ continue
def _AddSidechainFrameResidues(frame_residues, incomplete_sidechains,
keep_sidechains, res_list, rotamer_ids,
--
GitLab