IMP logo
IMP Reference Guide  2.17.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  """Keeps all structures inside a sphere."""
17 
18  def __init__(self, hierarchies, radius=10.0, resolution=10, weight=1.0,
19  center=None, label=None):
20  """Setup external barrier restraint.
21  @param hierarchies Can be one of the following inputs: IMP Hierarchy,
22  PMI System/State/Molecule/TempResidue, or a list/set of them
23  @param radius Size of external barrier
24  @param resolution Select which resolutions to act upon
25  @param weight Weight of restraint
26  @param center Center of the external barrier
27  (IMP.algebra.Vector3D object)
28  @param label A unique label to be used in outputs and
29  particle/restraint names.
30  """
31  hiers = IMP.pmi.tools.input_adaptor(hierarchies, resolution,
32  flatten=True)
33  model = hiers[0].get_model()
34  particles = [h.get_particle() for h in hiers]
35 
36  super(ExternalBarrier, self).__init__(model, label=label,
37  weight=weight)
38  self.radius = radius
39 
40  if center is None:
41  c3 = IMP.algebra.Vector3D(0, 0, 0)
42  elif type(center) is IMP.algebra.Vector3D:
43  c3 = center
44  else:
45  raise Exception(
46  "%s: @param center must be an IMP.algebra.Vector3D object" % (
47  self.name))
48 
49  ub3 = IMP.core.HarmonicUpperBound(radius, 10.0)
51  lsc = IMP.container.ListSingletonContainer(self.model)
52 
53  lsc.add(particles)
55  self.rs.add_restraint(r3)
56 
57 
58 class DistanceRestraint(IMP.pmi.restraints.RestraintBase):
59  """A simple distance restraint"""
60 
61  def __init__(self, root_hier, tuple_selection1, tuple_selection2,
62  distancemin=0, distancemax=100, resolution=1.0, kappa=1.0,
63  label=None, weight=1.):
64  """Setup distance restraint.
65  @param root_hier The hierarchy to select from
66  @param tuple_selection1 (resnum, resnum, molecule name, copy
67  number (=0))
68  @param tuple_selection2 (resnum, resnum, molecule name, copy
69  number (=0))
70  @param distancemin The minimum dist
71  @param distancemax The maximum dist
72  @param resolution For selecting particles
73  @param kappa The harmonic parameter
74  @param label A unique label to be used in outputs and
75  particle/restraint names
76  @param weight Weight of restraint
77  @note Pass the same resnum twice to each tuple_selection. Optionally
78  add a copy number.
79  """
80  ts1 = IMP.core.HarmonicUpperBound(distancemax, kappa)
81  ts2 = IMP.core.HarmonicLowerBound(distancemin, kappa)
82 
83  model = root_hier.get_model()
84  copy_num1 = 0
85  if len(tuple_selection1) > 3:
86  copy_num1 = tuple_selection1[3]
87  copy_num2 = 0
88  if len(tuple_selection2) > 3:
89  copy_num2 = tuple_selection2[3]
90 
91  sel1 = IMP.atom.Selection(root_hier,
92  resolution=resolution,
93  molecule=tuple_selection1[2],
94  residue_index=tuple_selection1[0],
95  copy_index=copy_num1)
96  particles1 = sel1.get_selected_particles()
97  sel2 = IMP.atom.Selection(root_hier,
98  resolution=resolution,
99  molecule=tuple_selection2[2],
100  residue_index=tuple_selection2[0],
101  copy_index=copy_num2)
102  particles2 = sel2.get_selected_particles()
103 
104  super(DistanceRestraint, self).__init__(model, label=label,
105  weight=weight)
106  print(self.name)
107 
108  print("Created distance restraint between "
109  "%s and %s" % (particles1[0].get_name(),
110  particles2[0].get_name()))
111 
112  if len(particles1) > 1 or len(particles2) > 1:
113  raise ValueError("more than one particle selected")
114 
115  self.rs.add_restraint(
116  IMP.core.DistanceRestraint(self.model, ts1,
117  particles1[0],
118  particles2[0]))
119  self.rs.add_restraint(
120  IMP.core.DistanceRestraint(self.model, ts2,
121  particles1[0],
122  particles2[0]))
123 
124 
126  """Restrain particles within (or outside) a cylinder.
127  The cylinder is aligned along the z-axis and with center x=y=0.
128  Optionally, one can restrain the cylindrical angle
129  """
130  import math
131 
132  def __init__(self, m, objects, resolution, radius, mintheta=None,
133  maxtheta=None, repulsive=False, label='None'):
134  '''
135  @param objects PMI2 objects to restrain
136  @param resolution the resolution you want the restraint to be applied
137  @param radius the radius of the cylinder
138  @param mintheta minimum cylindrical angle in degrees
139  @param maxtheta maximum cylindrical angle in degrees
140  @param repulsive If True, restrain the particles to be outside
141  of the cylinder instead of inside
142  @param label A unique label to be used in outputs and
143  particle/restraint names
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  '''Distance restraint with bistable potential
205  Authors: G. Bouvier, R. Pellarin. Pasteur Institute.
206  '''
207  import numpy as np
208  import math
209 
210  def __init__(self, m, p1, p2, dist1, dist2, sigma1, sigma2, weight1,
211  weight2):
212  '''
213  input two particles, the two equilibrium distances, their amplitudes,
214  and their weights (populations)
215  '''
216  IMP.Restraint.__init__(self, m, "BiStableDistanceRestraint %1%")
217  self.dist1 = dist1
218  self.dist2 = dist2
219 
220  self.sigma1 = sigma1
221  self.sigma2 = sigma2
222 
223  self.weight1 = weight1
224  self.weight2 = weight2
225 
226  if self.weight1+self.weight2 != 1:
227  raise ValueError("The sum of the weights must be one")
228 
229  self.d1 = IMP.core.XYZ(p1)
230  self.d2 = IMP.core.XYZ(p2)
231  self.particle_list = [p1, p2]
232 
233  def gaussian(self, x, mu, sig, w):
234  return (w*self.np.exp(-self.np.power(x - mu, 2.)
235  / (2 * self.np.power(sig, 2.))))
236 
237  def unprotected_evaluate(self, da):
238  dist = IMP.core.get_distance(self.d1, self.d2)
239  prob = self.gaussian(dist, self.dist1, self.sigma1, self.weight1) + \
240  self.gaussian(dist, self.dist2, self.sigma2, self.weight2)
241  return -self.math.log(prob)
242 
243  def do_get_inputs(self):
244  return self.particle_list
245 
246 
247 class DistanceToPointRestraint(IMP.pmi.restraints.RestraintBase):
248  """Anchor a particle to a specific coordinate."""
249 
250  def __init__(self, root_hier, tuple_selection,
251  anchor_point=IMP.algebra.Vector3D(0, 0, 0),
252  radius=10.0, kappa=10.0, resolution=1.0, weight=1.0,
253  label=None):
254  """Setup distance restraint.
255  @param root_hier The hierarchy to select from
256  @param tuple_selection (resnum, resnum, molecule name,
257  copy number (=0))
258  @param anchor_point Point to which to restrain particle
259  (IMP.algebra.Vector3D object)
260  @param radius Size of the tolerance length
261  @param kappa Strength of the harmonic restraint
262  @param resolution For selecting a particle
263  @param weight Weight of restraint
264  @param label A unique label to be used in outputs and
265  particle/restraint names
266  @note Pass the same resnum twice to each tuple_selection. Optionally
267  add a copy number
268  """
269  model = root_hier.get_model()
270  copy_num1 = 0
271  if len(tuple_selection) > 3:
272  copy_num1 = tuple_selection[3]
273 
274  sel1 = IMP.atom.Selection(root_hier,
275  resolution=resolution,
276  molecule=tuple_selection[2],
277  residue_index=tuple_selection[0],
278  copy_index=copy_num1)
279  ps = sel1.get_selected_particles()
280  if len(ps) > 1:
281  raise ValueError("More than one particle selected")
282 
283  super(DistanceToPointRestraint, self).__init__(model, label=label,
284  weight=weight)
285  self.radius = radius
286 
287  ub3 = IMP.core.HarmonicUpperBound(self.radius, kappa)
288  if anchor_point is None:
289  c3 = IMP.algebra.Vector3D(0, 0, 0)
290  elif isinstance(anchor_point, IMP.algebra.Vector3D):
291  c3 = anchor_point
292  else:
293  raise TypeError("anchor_point must be an algebra.Vector3D object")
294  ss3 = IMP.core.DistanceToSingletonScore(ub3, c3)
295 
296  lsc = IMP.container.ListSingletonContainer(self.model)
297  lsc.add(ps)
298 
299  r3 = IMP.container.SingletonsRestraint(ss3, lsc)
300  self.rs.add_restraint(r3)
301 
302  print("\n%s: Created distance_to_point_restraint between "
303  "%s and %s" % (self.name, ps[0].get_name(), c3))
304 
305 
306 class MembraneRestraint(IMP.pmi.restraints.RestraintBase):
307  """Restrain particles to be above, below, or inside a planar membrane.
308  The membrane is defined to lie on the xy plane with a given z
309  coordinate and thickness, and particles are restrained (by their
310  z coordinates) with a simple sigmoid score.
311  """
312 
313  def __init__(self, hier, objects_above=None, objects_inside=None,
314  objects_below=None, center=0.0, thickness=30.0,
315  softness=3.0, plateau=0.0000000001, resolution=1,
316  weight=1.0, label=None):
317  """Setup the restraint.
318 
319  @param objects_inside list or tuples of objects in membrane
320  (e.g. ['p1', (10, 30,'p2')])
321  @param objects_above list or tuples of objects above membrane
322  @param objects_below list or tuples of objects below membrane
323  @param thickness Thickness of the membrane along the z-axis
324  @param softness Softness of the limiter in the sigmoid function
325  @param plateau Parameter to set the probability (=1- plateau))
326  at the plateau phase of the sigmoid
327  @param weight Weight of restraint
328  @param label A unique label to be used in outputs and
329  particle/restraint names.
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  """Create an MRC density file to visualize the membrane."""
415  offset = 5.0 * self.thickness
416  apix = 3.0
417  resolution = 5.0
418 
419  # Create a density header of the requested size
421  IMP.algebra.Vector3D(-self.center - offset, -self.center - offset,
422  -self.center - offset),
423  IMP.algebra.Vector3D(self.center + offset, self.center + offset,
424  self.center + offset))
425  dheader = IMP.em.create_density_header(bbox, apix)
426  dheader.set_resolution(resolution)
427  dmap = IMP.em.SampledDensityMap(dheader)
428 
429  for vox in range(dmap.get_header().get_number_of_voxels()):
430  c = dmap.get_location_by_voxel(vox)
431  if self._is_membrane(c[2]) == 1:
432  dmap.set_value(c[0], c[1], c[2], 1.0)
433  else:
434  dmap.set_value(c[0], c[1], c[2], 0.0)
435 
436  IMP.em.write_map(dmap, file_out)
437 
438  def _is_membrane(self, z):
439  if ((z-self.center) < self.thickness/2.0 and
440  (z-self.center) >= -self.thickness/2.0):
441  return 1
442  else:
443  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
virtual double unprotected_evaluate(DerivativeAccumulator *da) const
Return the unweighted score for the restraint.
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:945
def __init__
input two 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.
Anchor a particle to a specific coordinate.
def create_membrane_density
Create an MRC density file to visualize the membrane.
def add_restraint_to_model
Add a PMI restraint to the model.
Definition: tools.py:91
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
Restrain particles to be above, below, or inside a planar membrane.
Keeps all structures inside a 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:48
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.
Distance restraint with bistable potential Authors: G.
Select hierarchy particles identified by the biological name.
Definition: Selection.h:70
virtual ModelObjectsTemp do_get_inputs() const =0
Restrain particles within (or outside) a cylinder.
A restraint is a term in an IMP ScoringFunction.
Definition: Restraint.h:53