IMP logo
IMP Reference Guide  2.15.0
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 
15 class ExternalBarrier(IMP.pmi.restraints.RestraintBase):
16 
17  """Restraint to keep all structures inside sphere."""
18 
19  def __init__(self, hierarchies, radius=10.0, resolution=10, weight=1.0,
20  center=None, label=None):
21  """Setup external barrier restraint.
22  @param radius Size of external barrier
23  @param hierarchies Can be one of the following inputs: IMP Hierarchy,
24  PMI System/State/Molecule/TempResidue, or a list/set of them
25  @param resolution Select which resolutions to act upon
26  @param weight Weight of restraint
27  @param center Center of the external barrier restraint
28  (IMP.algebra.Vector3D object)
29  @param label A unique label to be used in outputs and
30  particle/restraint names.
31  """
32  hiers = IMP.pmi.tools.input_adaptor(hierarchies, resolution,
33  flatten=True)
34  model = hiers[0].get_model()
35  particles = [h.get_particle() for h in hiers]
36 
37  super(ExternalBarrier, self).__init__(model, label=label,
38  weight=weight)
39  self.radius = radius
40 
41  if center is None:
42  c3 = IMP.algebra.Vector3D(0, 0, 0)
43  elif type(center) is IMP.algebra.Vector3D:
44  c3 = center
45  else:
46  raise Exception(
47  "%s: @param center must be an IMP.algebra.Vector3D object" % (
48  self.name))
49 
50  ub3 = IMP.core.HarmonicUpperBound(radius, 10.0)
52  lsc = IMP.container.ListSingletonContainer(self.model)
53 
54  lsc.add(particles)
56  self.rs.add_restraint(r3)
57 
58 
59 class DistanceRestraint(IMP.pmi.restraints.RestraintBase):
60 
61  """A simple distance restraint"""
62 
63  def __init__(self, root_hier, tuple_selection1, tuple_selection2,
64  distancemin=0, distancemax=100, resolution=1.0, kappa=1.0,
65  label=None, weight=1.):
66  """Setup distance restraint.
67  @param root_hier The hierarchy to select from
68  @param tuple_selection1 (resnum, resnum, molecule name, copy
69  number (=0))
70  @param tuple_selection2 (resnum, resnum, molecule name, copy
71  number (=0))
72  @param distancemin The minimum dist
73  @param distancemax The maximum dist
74  @param resolution For selecting particles
75  @param kappa The harmonic parameter
76  @param label A unique label to be used in outputs and
77  particle/restraint names
78  @param weight Weight of restraint
79  @note Pass the same resnum twice to each tuple_selection. Optionally
80  add a copy number (PMI2 only)
81  """
82  ts1 = IMP.core.HarmonicUpperBound(distancemax, kappa)
83  ts2 = IMP.core.HarmonicLowerBound(distancemin, kappa)
84 
85  model = root_hier.get_model()
86  copy_num1 = 0
87  if len(tuple_selection1) > 3:
88  copy_num1 = tuple_selection1[3]
89  copy_num2 = 0
90  if len(tuple_selection2) > 3:
91  copy_num2 = tuple_selection2[3]
92 
93  sel1 = IMP.atom.Selection(root_hier,
94  resolution=resolution,
95  molecule=tuple_selection1[2],
96  residue_index=tuple_selection1[0],
97  copy_index=copy_num1)
98  particles1 = sel1.get_selected_particles()
99  sel2 = IMP.atom.Selection(root_hier,
100  resolution=resolution,
101  molecule=tuple_selection2[2],
102  residue_index=tuple_selection2[0],
103  copy_index=copy_num2)
104  particles2 = sel2.get_selected_particles()
105 
106  super(DistanceRestraint, self).__init__(model, label=label,
107  weight=weight)
108  print(self.name)
109 
110  print("Created distance restraint between "
111  "%s and %s" % (particles1[0].get_name(),
112  particles2[0].get_name()))
113 
114  if len(particles1) > 1 or len(particles2) > 1:
115  raise ValueError("more than one particle selected")
116 
117  self.rs.add_restraint(
118  IMP.core.DistanceRestraint(self.model, ts1,
119  particles1[0],
120  particles2[0]))
121  self.rs.add_restraint(
122  IMP.core.DistanceRestraint(self.model, ts2,
123  particles1[0],
124  particles2[0]))
125 
126 
128  '''
129  PMI2 python restraint. Restrains particles within a
130  Cylinder aligned along the z-axis and
131  centered in x,y=0,0
132  Optionally, one can restrain the cylindrical angle
133  '''
134  import math
135 
136  def __init__(self, m, objects, resolution, radius, mintheta=None,
137  maxtheta=None, repulsive=False, label='None'):
138  '''
139  @param objects PMI2 objects
140  @param resolution the resolution you want the restraint to be applied
141  @param radius the radius of the cylinder
142  @param mintheta minimum cylindrical angle in degrees
143  @param maxtheta maximum cylindrical angle in degrees
144  '''
145  IMP.Restraint.__init__(self, m, "CylinderRestraint %1%")
146  self.radius = radius
147  self.softness = 3.0
148  self.softness_angle = 0.5
149  self.plateau = 1e-10
150  self.weight = 1.0
151  self.m = m
152  self.mintheta = mintheta
153  self.maxtheta = maxtheta
154  self.repulsive = repulsive
155  hierarchies = IMP.pmi.tools.input_adaptor(objects,
156  resolution,
157  flatten=True)
158  self.particles = [h.get_particle() for h in hierarchies]
159  self.label = label
160 
161  def get_probability(self, p):
162  xyz = IMP.core.XYZ(p)
163  r = self.math.sqrt(xyz.get_x()**2+xyz.get_y()**2)
164  argvalue = (r-self.radius) / self.softness
165  if self.repulsive:
166  argvalue = -argvalue
167  prob = (1.0 - self.plateau) / (1.0 + self.math.exp(-argvalue))
168  return prob
169 
170  def get_angle_probability(self, p):
171  xyz = IMP.core.XYZ(p)
172  angle = self.math.atan2(xyz.get_y(), xyz.get_x())*180.0/self.math.pi
173  anglediff = (angle - self.maxtheta + 180 + 360) % 360 - 180
174  argvalue1 = anglediff / self.softness_angle
175  anglediff = (angle - self.mintheta + 180 + 360) % 360 - 180
176  argvalue2 = -anglediff / self.softness_angle
177  prob = ((1.0-self.plateau)
178  / (1.0 + self.math.exp(-max(argvalue1, argvalue2))))
179  return prob
180 
181  def unprotected_evaluate(self, da):
182  s = 0.0
183  for p in self.particles:
184  s += -self.math.log(1.0-self.get_probability(p))
185  if self.mintheta is not None and self.maxtheta is not None:
186  s += -self.math.log(1.0-self.get_angle_probability(p))
187  return s
188 
189  def do_get_inputs(self):
190  return self.particles
191 
192  def add_to_model(self):
194 
195  def get_output(self):
196  output = {}
197  score = self.weight * self.unprotected_evaluate(None)
198  output["_TotalScore"] = str(score)
199  output["CylinderRestraint_" + self.label] = str(score)
200  return output
201 
202 
204  '''
205  a python restraint with bistable potential
206  Authors: G. Bouvier, R. Pellarin. Pasteur Institute.
207  '''
208  import numpy as np
209  import math
210 
211  def __init__(self, m, p1, p2, dist1, dist2, sigma1, sigma2, weight1,
212  weight2):
213  '''
214  input twp oarticles, the two equilibrium distances, their amplitudes,
215  and their weights (populations)
216  '''
217  IMP.Restraint.__init__(self, m, "BiStableDistanceRestraint %1%")
218  self.dist1 = dist1
219  self.dist2 = dist2
220 
221  self.sigma1 = sigma1
222  self.sigma2 = sigma2
223 
224  self.weight1 = weight1
225  self.weight2 = weight2
226 
227  if self.weight1+self.weight2 != 1:
228  raise ValueError("The sum of the weights must be one")
229 
230  self.d1 = IMP.core.XYZ(p1)
231  self.d2 = IMP.core.XYZ(p2)
232  self.particle_list = [p1, p2]
233 
234  def gaussian(self, x, mu, sig, w):
235  return (w*self.np.exp(-self.np.power(x - mu, 2.)
236  / (2 * self.np.power(sig, 2.))))
237 
238  def unprotected_evaluate(self, da):
239  dist = IMP.core.get_distance(self.d1, self.d2)
240  prob = self.gaussian(dist, self.dist1, self.sigma1, self.weight1) + \
241  self.gaussian(dist, self.dist2, self.sigma2, self.weight2)
242  return -self.math.log(prob)
243 
244  def do_get_inputs(self):
245  return self.particle_list
246 
247 
248 class DistanceToPointRestraint(IMP.pmi.restraints.RestraintBase):
249 
250  """Restraint for anchoring a particle to a specific coordinate."""
251 
252  def __init__(self, root_hier, tuple_selection,
253  anchor_point=IMP.algebra.Vector3D(0, 0, 0),
254  radius=10.0, kappa=10.0, resolution=1.0, weight=1.0,
255  label=None):
256  """Setup distance restraint.
257  @param root_hier The hierarchy to select from
258  @param tuple_selection (resnum, resnum, molecule name,
259  copy number (=0))
260  @param anchor_point Point to which to restrain particle
261  (IMP.algebra.Vector3D object)
262  @param radius Size of the tolerance length
263  @param kappa The harmonic parameter
264  @param resolution For selecting a particle
265  @param weight Weight of restraint
266  @param label A unique label to be used in outputs and
267  particle/restraint names
268  @note Pass the same resnum twice to each tuple_selection. Optionally
269  add a copy number (PMI2 only)
270  """
271  model = root_hier.get_model()
272  copy_num1 = 0
273  if len(tuple_selection) > 3:
274  copy_num1 = tuple_selection[3]
275 
276  sel1 = IMP.atom.Selection(root_hier,
277  resolution=resolution,
278  molecule=tuple_selection[2],
279  residue_index=tuple_selection[0],
280  copy_index=copy_num1)
281  ps = sel1.get_selected_particles()
282  if len(ps) > 1:
283  raise ValueError("More than one particle selected")
284 
285  super(DistanceToPointRestraint, self).__init__(model, label=label,
286  weight=weight)
287  self.radius = radius
288 
289  ub3 = IMP.core.HarmonicUpperBound(self.radius, kappa)
290  if anchor_point is None:
291  c3 = IMP.algebra.Vector3D(0, 0, 0)
292  elif isinstance(anchor_point, IMP.algebra.Vector3D):
293  c3 = anchor_point
294  else:
295  raise TypeError("anchor_point must be an algebra.Vector3D object")
296  ss3 = IMP.core.DistanceToSingletonScore(ub3, c3)
297 
298  lsc = IMP.container.ListSingletonContainer(self.model)
299  lsc.add(ps)
300 
301  r3 = IMP.container.SingletonsRestraint(ss3, lsc)
302  self.rs.add_restraint(r3)
303 
304  print("\n%s: Created distance_to_point_restraint between "
305  "%s and %s" % (self.name, ps[0].get_name(), c3))
306 
307 
308 class MembraneRestraint(IMP.pmi.restraints.RestraintBase):
309  def __init__(self, hier, objects_above=None, objects_inside=None,
310  objects_below=None, center=0.0, thickness=30.0,
311  softness=3.0, plateau=0.0000000001, resolution=1,
312  weight=1.0, label=None):
313  """ Setup Membrane restraint
314 
315  Simple sigmoid score calculated for particles above,
316 
317  @param objects_inside list or tuples of objects in membrane
318  (e.g. ['p1', (10, 30,'p2')])
319  @param objects_above list or tuples of objects above membrane
320  @param objects_below list or tuples of objects below membrane
321  @param thickness Thickness of the membrane along the z-axis
322  @param softness Softness of the limiter in the sigmoid function
323  @param plateau Parameter to set the probability (=1- plateau))
324  at the plateau phase of the sigmoid
325  @param weight Weight of restraint
326  @param label A unique label to be used in outputs and
327  particle/restraint names.
328  input a list of particles, the slope and theta of the sigmoid potential
329  theta is the cutoff distance for a protein-protein contact
330  """
331 
332  self.hier = hier
333  model = self.hier.get_model()
334 
335  super(MembraneRestraint, self).__init__(
336  model, name="MembraneRestraint", label=label, weight=weight)
337 
338  self.center = center
339  self.thickness = thickness
340  self.softness = softness
341  self.plateau = plateau
342  self.linear = 0.02
343  self.resolution = resolution
344 
345  # Create nuisance particle
346  p = IMP.Particle(model)
347  z_center = IMP.isd.Nuisance.setup_particle(p)
348  z_center.set_nuisance(self.center)
349 
350  # Setup restraint
351  mr = IMP.pmi.MembraneRestraint(model,
352  z_center.get_particle_index(),
353  self.thickness,
354  self.softness,
355  self.plateau,
356  self.linear)
357 
358  # Particles above
359  if objects_above:
360  for obj in objects_above:
361  if isinstance(obj, tuple):
362  self.particles_above = self._select_from_tuple(obj)
363 
364  elif isinstance(obj, str):
365  self.particles_above = self._select_from_string(obj)
366  mr.add_particles_above(self.particles_above)
367 
368  # Particles inside
369  if objects_inside:
370  for obj in objects_inside:
371  if isinstance(obj, tuple):
372  self.particles_inside = self._select_from_tuple(obj)
373 
374  elif isinstance(obj, str):
375  self.particles_inside = self._select_from_string(obj)
376  mr.add_particles_inside(self.particles_inside)
377 
378  # Particles below
379  if objects_below:
380  for obj in objects_below:
381  if isinstance(obj, tuple):
382  self.particles_below = self._select_from_tuple(obj)
383 
384  elif isinstance(obj, str):
385  self.particles_below = self._select_from_string(obj)
386  mr.add_particles_below(self.particles_below)
387 
388  self.rs.add_restraint(mr)
389 
390  def get_particles_above(self):
391  return self.particles_above
392 
393  def get_particles_inside(self):
394  return self.particles_inside
395 
396  def get_particles_below(self):
397  return self.particles_below
398 
399  def _select_from_tuple(self, obj):
400  particles = IMP.atom.Selection(
401  self.hier, molecule=obj[2],
402  residue_indexes=range(obj[0], obj[1]+1, 1),
403  resolution=self.resolution).get_selected_particles()
404 
405  return particles
406 
407  def _select_from_string(self, obj):
408  particles = IMP.atom.Selection(
409  self.hier, molecule=obj,
410  resolution=self.resolution).get_selected_particles()
411  return particles
412 
413  def create_membrane_density(self, file_out='membrane_localization.mrc'):
414 
415  '''
416  Just for visualization purposes.
417  Writes density of membrane localization
418  '''
419 
420  offset = 5.0 * self.thickness
421  apix = 3.0
422  resolution = 5.0
423 
424  # Create a density header of the requested size
426  IMP.algebra.Vector3D(-self.center - offset, -self.center - offset,
427  -self.center - offset),
428  IMP.algebra.Vector3D(self.center + offset, self.center + offset,
429  self.center + offset))
430  dheader = IMP.em.create_density_header(bbox, apix)
431  dheader.set_resolution(resolution)
432  dmap = IMP.em.SampledDensityMap(dheader)
433 
434  for vox in range(dmap.get_header().get_number_of_voxels()):
435  c = dmap.get_location_by_voxel(vox)
436  if self._is_membrane(c[2]) == 1:
437  dmap.set_value(c[0], c[1], c[2], 1.0)
438  else:
439  dmap.set_value(c[0], c[1], c[2], 0.0)
440 
441  IMP.em.write_map(dmap, file_out)
442 
443  def _is_membrane(self, z):
444  if ((z-self.center) < self.thickness/2.0 and
445  (z-self.center) >= -self.thickness/2.0):
446  return 1
447  else:
448  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)
Various classes to hold sets of particles.
Upper bound harmonic function (non-zero when feature > mean)
Miscellaneous utilities.
Definition: tools.py:1
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 input_adaptor
Adapt things for PMI (degrees of freedom, restraints, ...) Returns list of list of hierarchies...
Definition: tools.py:944
def __init__
input twp oarticles, 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:87
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.
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:421
Class to handle individual particles of a Model object.
Definition: Particle.h:41
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
virtual ModelObjectsTemp do_get_inputs() const =0
A restraint is a term in an IMP ScoringFunction.
Definition: Restraint.h:54