IMP logo
IMP Reference Guide  2.10.1
The Integrative Modeling Platform
restraints/basic.py
1 """@namespace IMP.pmi.restraints.basic
2 Some miscellaneous simple restraints.
3 """
4 
5 from __future__ import print_function
6 import IMP
7 import IMP.core
8 import IMP.algebra
9 import IMP.atom
10 import IMP.container
11 import IMP.pmi.tools
12 import IMP.pmi.restraints
13 
14 
16 
17  """Restraint to keep all structures inside sphere."""
18 
19  def __init__(self,
20  representation=None,
21  radius=10.0,
22  hierarchies=None,
23  resolution=10,
24  weight=1.0,
25  center=None,
26  label=None):
27  """Setup external barrier restraint.
28  @param representation DEPRECATED
29  @param radius Size of external barrier
30  @param hierarchies Can be one of the following inputs: IMP Hierarchy,
31  PMI System/State/Molecule/TempResidue, or a list/set of them
32  @param resolution Select which resolutions to act upon
33  @param weight Weight of restraint
34  @param center Center of the external barrier restraint
35  (IMP.algebra.Vector3D object)
36  @param label A unique label to be used in outputs and
37  particle/restraint names.
38  """
39  if representation:
40  model = representation.prot.get_model()
41  particles = IMP.pmi.tools.select(
42  representation,
43  resolution=resolution,
44  hierarchies=hierarchies)
45  elif hierarchies:
46  hiers = IMP.pmi.tools.input_adaptor(hierarchies, resolution,
47  flatten=True)
48  model = hiers[0].get_model()
49  particles = [h.get_particle() for h in hiers]
50  else:
51  raise Exception("%s: must pass representation or hierarchies" % (
52  self.name))
53 
54  super(ExternalBarrier, self).__init__(model, label=label,
55  weight=weight)
56  self.radius = radius
57 
58  if center is None:
59  c3 = IMP.algebra.Vector3D(0, 0, 0)
60  elif type(center) is IMP.algebra.Vector3D:
61  c3 = center
62  else:
63  raise Exception(
64  "%s: @param center must be an IMP.algebra.Vector3D object" % (
65  self.name))
66 
67  ub3 = IMP.core.HarmonicUpperBound(radius, 10.0)
69  lsc = IMP.container.ListSingletonContainer(self.model)
70 
71  lsc.add(particles)
73  self.rs.add_restraint(r3)
74 
75 
77 
78  """A simple distance restraint"""
79 
80  def __init__(self,
81  representation=None,
82  tuple_selection1=None,
83  tuple_selection2=None,
84  distancemin=0,
85  distancemax=100,
86  resolution=1.0,
87  kappa=1.0,
88  root_hier=None,
89  label=None,
90  weight=1.):
91  """Setup distance restraint.
92  @param representation DEPRECATED
93  @param tuple_selection1 (resnum, resnum, molecule name, copy
94  number (=0))
95  @param tuple_selection2 (resnum, resnum, molecule name, copy
96  number (=0))
97  @param distancemin The minimum dist
98  @param distancemax The maximum dist
99  @param resolution For selecting particles
100  @param kappa The harmonic parameter
101  @param root_hier The hierarchy to select from (use this instead of
102  representation)
103  @param label A unique label to be used in outputs and
104  particle/restraint names
105  @param weight Weight of restraint
106  @note Pass the same resnum twice to each tuple_selection. Optionally
107  add a copy number (PMI2 only)
108  """
109  if tuple_selection1 is None or tuple_selection2 is None:
110  raise Exception("You must pass tuple_selection1/2")
111  ts1 = IMP.core.HarmonicUpperBound(distancemax, kappa)
112  ts2 = IMP.core.HarmonicLowerBound(distancemin, kappa)
113 
114  if representation and not root_hier:
115  model = representation.prot.get_model()
116  particles1 = IMP.pmi.tools.select(representation,
117  resolution=resolution,
118  name=tuple_selection1[2],
119  residue=tuple_selection1[0])
120  particles2 = IMP.pmi.tools.select(representation,
121  resolution=resolution,
122  name=tuple_selection2[2],
123  residue=tuple_selection2[0])
124  elif root_hier and not representation:
125  model = root_hier.get_model()
126  copy_num1 = 0
127  if len(tuple_selection1) > 3:
128  copy_num1 = tuple_selection1[3]
129  copy_num2 = 0
130  if len(tuple_selection2) > 3:
131  copy_num2 = tuple_selection2[3]
132 
133  sel1 = IMP.atom.Selection(root_hier,
134  resolution=resolution,
135  molecule=tuple_selection1[2],
136  residue_index=tuple_selection1[0],
137  copy_index=copy_num1)
138  particles1 = sel1.get_selected_particles()
139  sel2 = IMP.atom.Selection(root_hier,
140  resolution=resolution,
141  molecule=tuple_selection2[2],
142  residue_index=tuple_selection2[0],
143  copy_index=copy_num2)
144  particles2 = sel2.get_selected_particles()
145  else:
146  raise Exception("Pass representation or root_hier, not both")
147 
148  super(DistanceRestraint, self).__init__(model, label=label,
149  weight=weight)
150  print(self.name)
151 
152  print("Created distance restraint between "
153  "%s and %s" % (particles1[0].get_name(),
154  particles2[0].get_name()))
155 
156  if len(particles1) > 1 or len(particles2) > 1:
157  raise ValueError("more than one particle selected")
158 
159  self.rs.add_restraint(
160  IMP.core.DistanceRestraint(self.model, ts1,
161  particles1[0],
162  particles2[0]))
163  self.rs.add_restraint(
164  IMP.core.DistanceRestraint(self.model, ts2,
165  particles1[0],
166  particles2[0]))
167 
168 
170  "If you use this class please let the PMI developers know.")
171 class TorqueRestraint(IMP.Restraint):
172  import math
173  def __init__(self, m, objects, resolution, angular_tolerance,label='None'):
174  IMP.Restraint.__init__(self, m, "TorqueRestraint %1%")
175  self.softness_angle = 0.5
176  self.plateau = 1e-10
177  self.weight = 1.0
178  self.m=m
179  hierarchies = IMP.pmi.tools.input_adaptor(objects,
180  resolution,
181  flatten=True)
182  self.particles = [h.get_particle() for h in hierarchies]
183  self.ds=[IMP.core.XYZ(p) for p in self.particles]
184  self.label=label
185  self.at=angular_tolerance
186 
187  def get_angle_probability(self,xyz,angle_center):
188  maxtheta=angle_center+self.at
189  mintheta=angle_center-self.at
190  angle=self.math.atan2(xyz.get_y(),xyz.get_x() )*180.0/self.math.pi
191  anglediff = (angle - maxtheta + 180 + 360) % 360 - 180
192  argvalue1=anglediff / self.softness_angle
193  anglediff = (angle - mintheta + 180 + 360) % 360 - 180
194  argvalue2=-anglediff / self.softness_angle
195  prob = (1.0-self.plateau) / (1.0 + self.math.exp(-max(argvalue1,argvalue2)))
196  return prob
197 
198  def unprotected_evaluate(self, da):
199  s=0.0
200  center=IMP.core.get_centroid(self.ds)
201  angle_center=self.math.atan2(center[1],center[0])*180.0/self.math.pi
202  for xyz in self.ds:
203  s+=-self.math.log(1.0-self.get_angle_probability(xyz,angle_center))
204  return s
205 
206  def do_get_inputs(self):
207  return self.particles
208 
209  def add_to_model(self):
211 
212  def get_output(self):
213  output = {}
214  score = self.weight * self.unprotected_evaluate(None)
215  output["_TotalScore"] = str(score)
216  output["TorqueRestraint_" + self.label] = str(score)
217  return output
218 
219 
220 
221 
223  '''
224  PMI2 python restraint. Restrains particles within a
225  Cylinder aligned along the z-axis and
226  centered in x,y=0,0
227  Optionally, one can restrain the cylindrical angle
228  '''
229  import math
230  def __init__(self, m, objects, resolution, radius,mintheta=None,
231  maxtheta=None,repulsive=False,label='None'):
232  '''
233  @param objects PMI2 objects
234  @param resolution the resolution you want the restraint to be applied
235  @param radius the radius of the cylinder
236  @param mintheta minimum cylindrical angle in degrees
237  @param maxtheta maximum cylindrical angle in degrees
238  '''
239  IMP.Restraint.__init__(self, m, "CylinderRestraint %1%")
240  self.radius=radius
241  self.softness = 3.0
242  self.softness_angle = 0.5
243  self.plateau = 1e-10
244  self.weight = 1.0
245  self.m=m
246  self.mintheta=mintheta
247  self.maxtheta=maxtheta
248  self.repulsive=repulsive
249  hierarchies = IMP.pmi.tools.input_adaptor(objects,
250  resolution,
251  flatten=True)
252  self.particles = [h.get_particle() for h in hierarchies]
253  self.label=label
254 
255  def get_probability(self,p):
256  xyz=IMP.core.XYZ(p)
257  r=self.math.sqrt(xyz.get_x()**2+xyz.get_y()**2)
258  argvalue=(r-self.radius) / self.softness
259  if self.repulsive: argvalue=-argvalue
260  prob = (1.0 - self.plateau) / (1.0 + self.math.exp(-argvalue))
261  return prob
262 
263  def get_angle_probability(self,p):
264  xyz=IMP.core.XYZ(p)
265  angle=self.math.atan2(xyz.get_y(),xyz.get_x() )*180.0/self.math.pi
266  anglediff = (angle - self.maxtheta + 180 + 360) % 360 - 180
267  argvalue1=anglediff / self.softness_angle
268  anglediff = (angle - self.mintheta + 180 + 360) % 360 - 180
269  argvalue2=-anglediff / self.softness_angle
270  prob = (1.0-self.plateau) / (1.0 + self.math.exp(-max(argvalue1,argvalue2)))
271  return prob
272 
273  def unprotected_evaluate(self, da):
274  s=0.0
275  for p in self.particles:
276  s+=-self.math.log(1.0-self.get_probability(p))
277  if self.mintheta is not None and self.maxtheta is not None:
278  s+=-self.math.log(1.0-self.get_angle_probability(p))
279  return s
280 
281  def do_get_inputs(self):
282  return self.particles
283 
284  def add_to_model(self):
286 
287  def get_output(self):
288  output = {}
289  score = self.weight * self.unprotected_evaluate(None)
290  output["_TotalScore"] = str(score)
291  output["CylinderRestraint_" + self.label] = str(score)
292  return output
293 
294 
295 
297  '''
298  a python restraint with bistable potential
299  Authors: G. Bouvier, R. Pellarin. Pasteur Institute.
300  '''
301  import numpy as np
302  import math
303 
304  def __init__(self,m,p1,p2,dist1,dist2,sigma1,sigma2,weight1,weight2):
305  '''
306  input twp particles, the two equilibrium distances, their amplitudes, and their weights (populations)
307  '''
308  IMP.Restraint.__init__(self, m, "BiStableDistanceRestraint %1%")
309  self.dist1=dist1
310  self.dist2=dist2
311 
312  self.sigma1=sigma1
313  self.sigma2=sigma2
314 
315  self.weight1=weight1
316  self.weight2=weight2
317 
318  if self.weight1+self.weight2 != 1:
319  raise ValueError("The sum of the weights must be one")
320 
321  self.d1=IMP.core.XYZ(p1)
322  self.d2=IMP.core.XYZ(p2)
323  self.particle_list=[p1,p2]
324 
325  def gaussian(self,x, mu, sig, w):
326  return w*self.np.exp(-self.np.power(x - mu, 2.) / (2 * self.np.power(sig, 2.)))
327 
328  def unprotected_evaluate(self,da):
329  dist=IMP.core.get_distance(self.d1,self.d2)
330  prob=self.gaussian(dist,self.dist1,self.sigma1,self.weight1)+\
331  self.gaussian(dist,self.dist2,self.sigma2,self.weight2)
332  return -self.math.log(prob)
333 
334  def do_get_inputs(self):
335  return self.particle_list
336 
337 
339 
340  """Restraint for anchoring a particle to a specific coordinate."""
341 
342  def __init__(self,
343  representation=None,
344  tuple_selection=None,
345  anchor_point=IMP.algebra.Vector3D(0, 0, 0),
346  radius=10.0,
347  kappa=10.0,
348  resolution=1.0,
349  weight=1.0,
350  root_hier=None,
351  label=None):
352  """Setup distance restraint.
353  @param representation DEPRECATED
354  @param tuple_selection (resnum, resnum, molecule name,
355  copy number (=0))
356  @param anchor_point Point to which to restrain particle
357  (IMP.algebra.Vector3D object)
358  @param radius Size of the tolerance length
359  @param kappa The harmonic parameter
360  @param resolution For selecting a particle
361  @param weight Weight of restraint
362  @param root_hier The hierarchy to select from (use this instead of
363  representation)
364  @param label A unique label to be used in outputs and
365  particle/restraint names
366  @note Pass the same resnum twice to each tuple_selection. Optionally
367  add a copy number (PMI2 only)
368  """
369  if tuple_selection is None:
370  raise ValueError("You must pass a tuple_selection")
371 
372  if representation and not root_hier:
373  model = representation.prot.get_model()
374  ps = IMP.pmi.tools.select(representation,
375  resolution=resolution,
376  name=tuple_selection[2],
377  residue=tuple_selection[0])
378  elif root_hier and not representation:
379  model = root_hier.get_model()
380  copy_num1 = 0
381  if len(tuple_selection) > 3:
382  copy_num1 = tuple_selection[3]
383 
384  sel1 = IMP.atom.Selection(root_hier,
385  resolution=resolution,
386  molecule=tuple_selection[2],
387  residue_index=tuple_selection[0],
388  copy_index=copy_num1)
389  ps = sel1.get_selected_particles()
390  else:
391  raise ValueError("Pass representation or root_hier, not both")
392  if len(ps) > 1:
393  raise ValueError("More than one particle selected")
394 
395  super(DistanceToPointRestraint, self).__init__(model, label=label,
396  weight=weight)
397  self.radius = radius
398 
399  ub3 = IMP.core.HarmonicUpperBound(self.radius, kappa)
400  if anchor_point is None:
401  c3 = IMP.algebra.Vector3D(0, 0, 0)
402  elif isinstance(anchor_point, IMP.algebra.Vector3D):
403  c3 = anchor_point
404  else:
405  raise TypeError("anchor_point must be an algebra.Vector3D object")
406  ss3 = IMP.core.DistanceToSingletonScore(ub3, c3)
407 
408  lsc = IMP.container.ListSingletonContainer(self.model)
409  lsc.add(ps)
410 
411  r3 = IMP.container.SingletonsRestraint(ss3, lsc)
412  self.rs.add_restraint(r3)
413 
414  print("\n%s: Created distance_to_point_restraint between "
415  "%s and %s" % (self.name, ps[0].get_name(), c3))
416 
418  def __init__(self,
419  hier,
420  objects_above=None,
421  objects_inside=None,
422  objects_below=None,
423  center = 0.0,
424  thickness=30.0,
425  softness=3.0,
426  plateau=0.0000000001,
427  resolution=1,
428  weight = 1.0,
429  label = None):
430  """ Setup Membrane restraint
431 
432  Simple sigmoid score calculated for particles above,
433 
434  @param objects_inside list or tuples of objects in membrane (e.g. ['p1', (10, 30,'p2')])
435  @param objects_above list or tuples of objects above membrane
436  @param objects_below list or tuples of objects below membrane
437  @param thickness Thickness of the membrane along the z-axis
438  @param softness Softness of the limiter in the sigmoid function
439  @param plateau Parameter to set the probability (=1- plateau)) at the plateau
440  phase of the sigmoid
441  @param weight Weight of restraint
442  @param label A unique label to be used in outputs and
443  particle/restraint names.
444  input a list of particles, the slope and theta of the sigmoid potential
445  theta is the cutoff distance for a protein-protein contact
446  """
447 
448  self.hier = hier
449  model = self.hier.get_model()
450 
451  rname = "MembraneRestraint"
452  super(MembraneRestraint, self).__init__(
453  model, name="MembraneRestraint", label=label, weight=weight)
454 
455  self.center = center
456  self.thickness = thickness
457  self.softness = softness
458  self.plateau = plateau
459  self.linear = 0.02
460  self.resolution = resolution
461 
462  # Create nuisance particle
463  p = IMP.Particle(model)
464  z_center = IMP.isd.Nuisance.setup_particle(p)
465  z_center.set_nuisance(self.center)
466 
467  # Setup restraint
468  mr = IMP.pmi.MembraneRestraint(model,
469  z_center.get_particle_index(),
470  self.thickness,
471  self.softness,
472  self.plateau,
473  self.linear)
474 
475  # Particles above
476  if objects_above:
477  for obj in objects_above:
478  if isinstance(obj, tuple):
479  self.particles_above = self._select_from_tuple(obj)
480 
481  elif isinstance(obj, str):
482  self.particles_above = self._select_from_string(obj)
483  mr.add_particles_above(self.particles_above)
484 
485  # Particles inside
486  if objects_inside:
487  for obj in objects_inside:
488  if isinstance(obj, tuple):
489  self.particles_inside = self._select_from_tuple(obj)
490 
491  elif isinstance(obj, str):
492  self.particles_inside = self._select_from_string(obj)
493  mr.add_particles_inside(self.particles_inside)
494 
495 
496  # Particles below
497  if objects_below:
498  for obj in objects_below:
499  if isinstance(obj, tuple):
500  self.particles_below = self._select_from_tuple(obj)
501 
502  elif isinstance(obj, str):
503  self.particles_below = self._select_from_string(obj)
504  mr.add_particles_below(self.particles_below)
505 
506 
507  self.rs.add_restraint(mr)
508 
509  def get_particles_above(self):
510  return self.particles_above
511 
512  def get_particles_inside(self):
513  return self.particles_inside
514 
515  def get_particles_below(self):
516  return self.particles_below
517 
518  def _select_from_tuple(self, obj):
519  particles = IMP.atom.Selection(self.hier,
520  molecule = obj[2],
521  residue_indexes = range(obj[0], obj[1]+1, 1),
522  resolution = self.resolution).get_selected_particles()
523 
524  return particles
525 
526  def _select_from_string(self, obj):
527  particles = IMP.atom.Selection(self.hier,
528  molecule = obj,
529  resolution = self.resolution).get_selected_particles()
530  return particles
531 
532  def create_membrane_density(self, file_out='membrane_localization.mrc'):
533 
534  '''
535  Just for visualization purposes.
536  Writes density of membrane localization
537  '''
538 
539  offset = 5.0*self.thickness
540  apix = 3.0
541  resolution = 5.0
542 
543  # Create a density header of the requested size
545  IMP.algebra.Vector3D(-self.center - offset, -self.center - offset, -self.center - offset,),
546  IMP.algebra.Vector3D(self.center + offset, self.center + offset, self.center + offset))
547  dheader = IMP.em.create_density_header(bbox, apix)
548  dheader.set_resolution(resolution)
549  dmap = IMP.em.SampledDensityMap(dheader)
550 
551  for vox in range(dmap.get_header().get_number_of_voxels()):
552  c = dmap.get_location_by_voxel(vox)
553  if self._is_membrane(c[2])==1:
554  dmap.set_value(c[0], c[1], c[2], 1.0)
555  else:
556  dmap.set_value(c[0], c[1], c[2], 0.0)
557 
558  IMP.em.write_map(dmap, file_out)
559 
560  def _is_membrane(self, z):
561  if (z-self.center) < self.thickness/2.0 and (z-self.center) >= -self.thickness/2.0 :
562  return 1
563  else:
564  return 0
Applies a SingletonScore to each Singleton in a list.
static Nuisance setup_particle(Model *m, ParticleIndex pi)
Definition: Nuisance.h:31
def __init__
Setup distance restraint.
Lower bound harmonic function (non-zero when feature < mean)
def add_to_model
Add the restraint to the model.
Various classes to hold sets of particles.
Upper bound harmonic function (non-zero when feature > mean)
Miscellaneous utilities.
Definition: tools.py:1
algebra::Vector3D get_centroid(const XYZs &ps)
Get the centroid.
Distance restraint between two particles.
double get_distance(XYZR a, XYZR b)
Compute the sphere distance between a and b.
Definition: XYZR.h:89
Classes to handle different kinds of restraints.
def __init__
input twp particles, the two equilibrium distances, their amplitudes, and their weights (populations)...
Class for sampling a density map from particles.
def __init__
Setup external barrier restraint.
Store a list of ParticleIndexes.
Restraint for anchoring a particle to a specific coordinate.
def add_restraint_to_model
Add a PMI restraint to the model.
Definition: tools.py:69
DensityHeader create_density_header(const algebra::BoundingBoxD< 3 > &bb, float spacing)
Create a header from a bounding box in 3D.
Apply a function to the distance to a fixed point.
A decorator for a particle with x,y,z coordinates.
Definition: XYZ.h:30
Restraint to keep all structures inside sphere.
def deprecated_object
Python decorator to mark a class as deprecated.
Definition: __init__.py:9687
Basic functionality that is expected to be used by a wide variety of IMP users.
General purpose algebraic and geometric methods that are expected to be used by a wide variety of IMP...
The general base class for IMP exceptions.
Definition: exception.h:49
VectorD< 3 > Vector3D
Definition: VectorD.h:395
Class to handle individual particles of a Model object.
Definition: Particle.h:41
def select
this function uses representation=SimplifiedModel it returns the corresponding selected particles rep...
Definition: tools.py:726
Functionality for loading, creating, manipulating and scoring atomic structures.
a python restraint with bistable potential Authors: G.
Select hierarchy particles identified by the biological name.
Definition: Selection.h:66
def get_output
Get outputs to write to stat files.
virtual ModelObjectsTemp do_get_inputs() const =0
Base class for PMI restraints, which wrap IMP.Restraint(s).
A restraint is a term in an IMP ScoringFunction.
Definition: Restraint.h:54