OpenStructure
helix_kinks.py
Go to the documentation of this file.
1 """
2 Functions to calculate helix kinks: bend, face shift and wobble angles
3 
4 Author: Niklaus Johner
5 """
6 
7 import os
8 import ost
9 
10 def __FindProline(sele,proline):
11  if not sele.IsValid():
12  print('selection is not valid')
13  raise RuntimeError
14  if proline==False:
15  proline=sele.Select('rname=PRO')
16  if not proline.GetResidueCount()==1:
17  print(proline.GetResidueCount(),'prolines in the selection. One proline is needed')
18  raise RuntimeError
19  proline=proline.residues[0]
20  proline_ca=proline.FindAtom('CA')
21  if not proline_ca.IsValid():
22  print('proline has no CA atom')
23  raise RuntimeError
24  return (proline,proline_ca)
25  proline.GetNumber().num
26 
27 def __SelectPreAndPostProline(sele,proline_num):
28  pre_proline=sele.Select('rnum<'+str(proline_num))
29  post_proline=sele.Select('rnum>'+str(proline_num))
30  print('pre-proline residues')
31  for res in pre_proline.residues:
32  print(res)
33  print('post-proline residues')
34  for res in post_proline.residues:
35  print(res)
36  if pre_proline.GetResidueCount()<4 or post_proline.GetResidueCount()<4:
37  print('pre and post proline helices should be at least 4 residues long, 7 for better stability')
38  raise RuntimeError
39  return (pre_proline,post_proline)
40 
41 def __FindCa3AndCa4(sele,proline_ca,proline_num):
42  ca_3=sele.FindAtom(proline_ca.GetHandle().GetChain().GetName(),proline_num-3,'CA')
43  ca_4=sele.FindAtom(proline_ca.GetHandle().GetChain().GetName(),proline_num-4,'CA')
44  if not (ca_3.IsValid() and ca_4.IsValid()):
45  print('CA not found in (i-4) or (i-3) residue')
46  raise RuntimeError
47  return (ca_3,ca_4)
48 
49 
50 def __CalculateBendAngle(pre_proline_axis,post_proline_axis):
51  return ost.geom.Angle(pre_proline_axis,post_proline_axis)
52 
53 def __CalculateWobbleAngle(pre_proline_axis,post_proline_axis,post_proline_centers,proline_pos):
54  p1=proline_pos-post_proline_centers
55  n1=p1-ost.geom.Dot(p1,post_proline_axis)*post_proline_axis
56  p2=-pre_proline_axis
57  n2=p2-ost.geom.Dot(p2,post_proline_axis)*post_proline_axis
58  sign=ost.geom.Dot(ost.geom.Cross(pre_proline_axis,n2),n2-n1)
59  sign=sign/abs(sign)
60  return sign*ost.geom.Angle(n1,n2)
61 
62 def __CalculateFaceShift(pre_proline_axis,post_proline_axis,pre_proline_centers,post_proline_centers,proline_pos,ca3_pos,ca4_pos,bend_angle):
63  p1=proline_pos-post_proline_centers
64  n1=p1-ost.geom.Dot(p1,post_proline_axis)*post_proline_axis
65  p2=(ca3_pos+ca4_pos)/2.0-pre_proline_centers
66  n2=p2-ost.geom.Dot(p2,pre_proline_axis)*pre_proline_axis
67  #Here we want to apply a rotation to n1, to align the post-proline axis on the pre-proline axis
68  R=ost.geom.AxisRotation(ost.geom.Cross(post_proline_axis,pre_proline_axis),bend_angle)
69  n1=R*n1
70  #We also need to determine the sign of the angle
71  sign=ost.geom.Dot(ost.geom.Cross(pre_proline_axis,n2),n2-n1)
72  sign=sign/abs(sign)
73  return sign*ost.geom.Angle(n1,n2)
74 
75 
76 def AnalyzeHelixKink(t, sele, proline=False):
77  """
78  This function calculates the bend, wobble and face-shift angles in an alpha-
79  helix over a trajectory. The determination is more stable if there are at
80  least 4 residues on each side (8 is even better) of the proline around which
81  the helix is kinked. The selection should contain all residues in the correct
82  order and with no gaps and no missing C-alphas.
83 
84  :param t: The trajectory to be analyzed
85  :type t: :class:`~ost.mol.CoordGroup`
86  :param sele: A selection containing the alpha helix to be analyzed
87  :type sele: :class:`~ost.mol.EntityView`
88  :param proline: A selection containing only the proline (or another residue)
89  around which the helix is kinked. If False, the proline will
90  be searched for automatically
91  :type proline: :class:`ost.mol.EntityView`
92 
93  :return: A tuple (bend_angle, face_shift, wobble_angle).
94  :rtype: (FloatList, FLoatList, FloatList)
95  """
96  n_frames=t.GetFrameCount()
97  (proline,proline_ca)=__FindProline(sele,proline)
98  proline_num=proline.GetNumber().num
99  (pre_proline,post_proline)=__SelectPreAndPostProline(sele,proline_num)
100  (ca_3,ca_4)=__FindCa3AndCa4(sele,proline_ca,proline_num)
101  #Here we extract the necessary information from the trajectory
102  pre_proline_axis=ost.geom.Vec3List()
103  post_proline_axis=ost.geom.Vec3List()
104  pre_proline_centers=ost.geom.Vec3List()
105  post_proline_centers=ost.geom.Vec3List()
106  ost.mol.alg.AnalyzeAlphaHelixAxis(t,pre_proline,pre_proline_axis,pre_proline_centers)
107  ost.mol.alg.AnalyzeAlphaHelixAxis(t,post_proline,post_proline_axis,post_proline_centers)
108  proline_pos=ost.mol.alg.AnalyzeAtomPos(t,proline_ca.GetHandle())
109  ca3_pos=ost.mol.alg.AnalyzeAtomPos(t,ca_3.GetHandle())
110  ca4_pos=ost.mol.alg.AnalyzeAtomPos(t,ca_4.GetHandle())
111  #Now we calculate the bend angle
112  bend_angle=ost.FloatList()
113  face_shift=ost.FloatList()
114  wobble_angle=ost.FloatList()
115  for i in range(n_frames):
116  bend_angle.append(__CalculateBendAngle(pre_proline_axis[i],post_proline_axis[i]))
117  face_shift.append(__CalculateFaceShift(pre_proline_axis[i],post_proline_axis[i],pre_proline_centers[i],post_proline_centers[i],proline_pos[i],ca3_pos[i],ca4_pos[i],bend_angle[i]))
118  wobble_angle.append(__CalculateWobbleAngle(pre_proline_axis[i],post_proline_axis[i],post_proline_centers[i],proline_pos[i]))
119  return (bend_angle,face_shift,wobble_angle)
120 
121 
122 def CalculateHelixKink(sele, proline=False):
123  """
124  This function calculates the bend, wobble and face-shift angles in an alpha-
125  helix of an EntityView. The determination is more stable if there are at least
126  4 residues on each side (8 is even better) of the proline around which the
127  helix is kinked. The selection should contain all residues in the correct
128  order and with no gaps and no missing C-alphas.
129 
130  :param sele: A selection containing the alpha helix to be analyzed
131  :type sele: :class:`~ost.mol.EntityView`
132  :param proline: A selection containing only the proline (or another residue)
133  around which the helix is kinked. If False, the proline will
134  be searched for automatically
135  :type proline: :class:`ost.mol.EntityView`
136 
137  :return: A tuple (bend_angle, face_shift, wobble_angle).
138  :rtype: (float, float, float)
139  """
140  (proline,proline_ca)=__FindProline(sele,proline)
141  proline_num=proline.GetNumber().num
142  (pre_proline,post_proline)=__SelectPreAndPostProline(sele,proline_num)
143  (ca_3,ca_4)=__FindCa3AndCa4(sele,proline_ca,proline_num)
144  #Here we extract the necessary information from the structure
145  pre_proline_axis=ost.mol.alg.structure_analysis.CalculateHelixAxis(pre_proline)
146  post_proline_axis=ost.mol.alg.structure_analysis.CalculateHelixAxis(post_proline)
147  prepa=pre_proline_axis.GetDirection()
148  prepc=pre_proline_axis.GetOrigin()
149  postpa=post_proline_axis.GetDirection()
150  postpc=post_proline_axis.GetOrigin()
151  #calculate the angles
152  bend=__CalculateBendAngle(prepa,postpa)
153  wobble=__CalculateWobbleAngle(prepa,postpa,postpc,proline_ca.pos)
154  shift=__CalculateFaceShift(prepa,postpa,prepc,postpc,proline_ca.pos,ca_3.pos,ca_4.pos,bend)
155  return (bend,shift,wobble)
156 
157 
158 
159 
def CalculateHelixKink(sele, proline=False)
Definition: helix_kinks.py:122
def AnalyzeHelixKink(t, sele, proline=False)
Definition: helix_kinks.py:76
void DLLEXPORT_OST_MOL_ALG AnalyzeAlphaHelixAxis(const CoordGroupHandle &traj, const EntityView &prot_seg, geom::Vec3List &directions, geom::Vec3List &centers, unsigned int stride=1)
geom::Vec3List DLLEXPORT_OST_MOL_ALG AnalyzeAtomPos(const CoordGroupHandle &traj, const AtomHandle &a1, unsigned int stride=1)