OpenStructure
 All Data Structures Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Groups Pages
ialign.py
Go to the documentation of this file.
1 #------------------------------------------------------------------------------
2 # This file is part of the OpenStructure project <www.openstructure.org>
3 #
4 # Copyright (C) 2008-2009 by the OpenStructure authors
5 #
6 # This library is free software; you can redistribute it and/or modify it under
7 # the terms of the GNU Lesser General Public License as published by the Free
8 # Software Foundation; either version 3.0 of the License, or (at your option)
9 # any later version.
10 # This library is distributed in the hope that it will be useful, but WITHOUT
11 # ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
12 # FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License for more
13 # details.
14 #
15 # You should have received a copy of the GNU Lesser General Public License
16 # along with this library; if not, write to the Free Software Foundation, Inc.,
17 # 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA
18 #------------------------------------------------------------------------------
19 """
20 Wrappers for the iAlign utility.
21 
22 References:
23 
24 Mu Gao and Jeffrey Skolnick, 2010. iAlign: a method for the structural comparison
25 of protein-protein interfaces. Bioinformatics, 26(18):2259-65.
26 Mu Gao and Jeffrey Skolnick, 2010. Structural space of protein-protein interfaces
27 is degenerate, close to complete, and highly connected. PNAS 107(52):22517-22.
28 
29 Authors: Pascal Benkert, Marco Biasini, Martino Bertoni
30 """
31 
32 import subprocess, os, tempfile, platform
33 from ost import settings, io, geom, seq
34 
35 def _SetupFiles(models):
36  # create temporary directory
37  tmp_dir_name=tempfile.mkdtemp()
38  dia = 'PDB'
39  for index, model in enumerate(models):
40  for chain in model.chains:
41  if len(chain.name) > 1:
42  dia = 'CHARMM'
43  break;
44  for res in chain.residues:
45  if len(res.name) > 3:
46  dia = 'CHARMM'
47  break;
48  io.SavePDB(model, os.path.join(tmp_dir_name, 'model%02d.pdb' % (index+1)), dialect=dia)
49  return tmp_dir_name
50 
51 def _CleanupFiles(dir_name):
52  import shutil
53  shutil.rmtree(dir_name)
54 
56  """
57  Holds the result of running iAlign
58 
59  .. attribute:: rmsd
60 
61  The RMSD of the common Calpha atoms of both structures
62 
63  .. attribute:: transform
64 
65  The transform that superposes the model onto the reference structure.
66 
67  :type: :class:`~ost.geom.Mat4`
68 
69  .. attribute:: alignment
70 
71  The alignment of the structures, that is the pairing of Calphas of both
72  structures. Since the programs only read ATOM records, residues consisting
73  of HETATMs (MSE) are not included in the alignment.
74 
75  :type: :class:`~ost.seq.AlignmentHandle`
76 
77  .. attribute:: IS_score
78 
79  The IS-score of the structural superposition
80 
81  .. attribute:: aligned_residues
82 
83  The total number of aligned residues
84 
85  .. attribute:: aligned_contacts
86 
87  The total number of aligned contacts
88 
89  """
90  def __init__(self, rmsd, transform, alignment, is_score,
91  aligned_residues, aligned_contacts):
92 
93  self.rmsd=rmsd
94  self.transform=transform
95  self.alignment=alignment
96  self.is_score=is_score
97  self.aligned_residues=aligned_residues
98  self.aligned_contacts=aligned_contacts
99 
100 def _ParseiAlign(lines):
101  info_line=lines[18].split(',')
102  is_score=float(info_line[0].split('=')[1].strip())
103  aln_residues=int(lines[19].split('=')[1].strip())
104  aln_contacts=int(lines[20].split('=')[1].strip())
105  info_line=lines[21].split(',')
106  rmsd=float(info_line[0].split('=')[1].strip())
107 
108  tf1=[float(i.strip()) for i in lines[25][1:].split()]
109  tf2=[float(i.strip()) for i in lines[26][1:].split()]
110  tf3=[float(i.strip()) for i in lines[27][1:].split()]
111  rot=geom.Mat3(tf1[2], tf1[3], tf1[4], tf2[2], tf2[3],
112  tf2[4], tf3[2], tf3[3], tf3[4])
113  tf=geom.Mat4(rot)
114  tf.PasteTranslation(geom.Vec3(tf1[1], tf2[1], tf3[1]))
115  seq1 = seq.CreateSequence("1",lines[32].strip())
116  seq2 = seq.CreateSequence("2",lines[34].strip())
117  alignment = seq.CreateAlignment()
118  alignment.AddSequence(seq2)
119  alignment.AddSequence(seq1)
120  return iAlignResult(rmsd, tf, alignment, is_score, aln_residues, aln_contacts)
121 
122 def _RuniAlign(ialign, tmp_dir, options={}):
123  opts = {'a' : 1, # concise output
124  'w' : tmp_dir
125  }
126  opts.update(options)
127  cmd_opts = []
128  for k, v in opts.iteritems():
129  if type(v) == type(True):
130  if v == True:
131  cmd_opts.append('-%s' % str(k))
132  else:
133  cmd_opts.append('-%s %s' % (str(k), str(v)))
134  cmd_opts = ' '.join(cmd_opts)
135  model1_filename=os.path.join(tmp_dir, 'model01.pdb')
136  model2_filename=os.path.join(tmp_dir, 'model02.pdb')
137  if platform.system() == "Windows":
138  ialign_path=settings.Locate('ialign.pl', explicit_file_name=ialign)
139  command="\"%s\" %s %s %s" % (os.path.normpath(ialign_path), model1_filename, model2_filename, cmd_opts)
140  else:
141  ialign_path=settings.Locate('ialign.pl', explicit_file_name=ialign)
142  command="\"%s\" \"%s\" \"%s\" %s" % (ialign_path, model1_filename, model2_filename, cmd_opts)
143  ps=subprocess.Popen(command, shell=True, stdout=subprocess.PIPE)
144  ps.wait()
145  lines=ps.stdout.readlines()
146  if (len(lines))<22:
147  _CleanupFiles(tmp_dir)
148  #for l in lines:
149  # print l
150  raise RuntimeError("iAlign superposition failed")
151  return _ParseiAlign(lines)
152 
153 
154 def iAlign(model1, model2, ialign=None):
155  """
156  Compare protein-protein interfaces of the structures of two pairs of
157  protein complexes and suporpose them.
158 
159 
160  :param model1: The model structure. If the superposition is successful, will
161  be superposed onto the reference structure
162  :type model1: :class:`~ost.mol.EntityView` or :class:`~ost.mol.EntityHandle`
163  :param model2: The reference structure
164  :type model2: :class:`~ost.mol.EntityView` or :class:`~ost.mol.EntityHandle`
165  :param ialign: If not None, the path to the ialign executable.
166  :returns: The result of the tmscore superposition
167  :rtype: :class:`iAlignResult`
168 
169  :raises: :class:`~ost.settings.FileNotFound` if tmalign could not be located.
170  :raises: :class:`RuntimeError` if the superposition failed
171  """
172  tmp_dir_name=_SetupFiles((model1, model2))
173  result=_RuniAlign(ialign, tmp_dir_name)
174  model1.handle.EditXCS().ApplyTransform(result.transform)
175  _CleanupFiles(tmp_dir_name)
176  return result
Three dimensional vector class, using Real precision.
Definition: vec3.hh:42