IMP logo
IMP Reference Guide  2.10.1
The Integrative Modeling Platform
pmi/restraints/proteomics.py
1 """@namespace IMP.pmi.restraints.proteomics
2 Restraints for handling various kinds of proteomics data.
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
12 import IMP.pmi.tools
13 import IMP.pmi.output
14 import numpy
15 import math
16 import math
17 import sys
18 
19 
20 class ConnectivityRestraint(object):
21 
22  '''
23  generate a connectivity restraint between domains
24  setting up the restraint
25  example:
26  cr=restraints.ConnectivityRestraint(simo,["CCC",(1,100,"TTT"),(100,150,"AAA")])
27  cr.add_to_model()
28  cr.set_label("CR1")
29  Multistate support =No
30  Selection type=selection tuple
31  Resolution=Yes
32  '''
33 
34  def __init__(
35  self,
36  representation,
37  selection_tuples,
38  kappa=10.0,
39  resolution=None,
40  label="None"):
41 
42  self.weight = 1.0
43  self.kappa = kappa
44  self.label = label
45  if self.label == "None":
46  self.label = str(selection_tuples)
47  self.m = representation.prot.get_model()
48  self.rs = IMP.RestraintSet(self.m, label)
49 
50  sels = []
51 
52  for s in selection_tuples:
53  particles = IMP.pmi.tools.select_by_tuple(representation, s,
54  resolution=resolution, name_is_ambiguous=True)
55  sel = IMP.atom.Selection(particles)
56  sels.append(sel)
57 
59  sels,
60  self.kappa,
61  self.label)
62  self.rs.add_restraint(cr)
63 
64  def set_label(self, label):
65  self.label = label
66  self.rs.set_name(label)
67  for r in self.rs.get_restraints():
68  r.set_name(label)
69 
70  def add_to_model(self):
72 
73  def get_restraint(self):
74  return self.rs
75 
76  def get_restraints(self):
77  rlist = []
78  for r in self.rs.get_restraints():
79  rlist.append(IMP.core.PairRestraint.get_from(r))
80  return rlist
81 
82  def set_weight(self, weight):
83  self.weight = weight
84  self.rs.set_weight(weight)
85 
86  def get_output(self):
87  output = {}
88  score = self.weight * self.rs.unprotected_evaluate(None)
89  output["_TotalScore"] = str(score)
90  output["ConnectivityRestraint_" + self.label] = str(score)
91  return output
92 
93 #
94 
95 
96 class CompositeRestraint(object):
97 
98  '''
99  handleparticles is a selection tuple
100  compositeparticles is a list of selection tuples
101  '''
102 
103  def __init__(
104  self,
105  representation,
106  handleparticles_tuples,
107  compositeparticles_tuple_list,
108  cut_off=5.0,
109  lam=1.0,
110  plateau=0.0,
111  resolution=None,
112  label="None"):
113 
114  # composite particles: all particles beside the handle
115  self.label = label
116  self.m = representation.prot.get_model()
117  self.rs = IMP.RestraintSet(self.m, 'cr')
118 
119  self.handleparticles = []
120  for s in handleparticles_tuples:
121  self.handleparticles += IMP.pmi.tools.select_by_tuple(
122  representation, s,
123  resolution=resolution, name_is_ambiguous=True)
124  self.compositeparticles = []
125  compositeparticle_list = []
126  for list in compositeparticles_tuple_list:
127  tmplist = []
128  for s in list:
129  tmplist += IMP.pmi.tools.select_by_tuple(
130  representation, s,
131  resolution=resolution, name_is_ambiguous=True)
132  compositeparticle_list.append(tmplist)
133  self.compositeparticles += tmplist
134 
136  self.m,
137  self.handleparticles,
138  cut_off,
139  lam,
140  True,
141  plateau)
142 
143  for ps in compositeparticle_list:
144  # composite particles is a list of list of particles
145  ln.add_composite_particle(ps)
146 
147  self.rs.add_restraint(ln)
148 
149  def set_label(self, label):
150  self.label = label
151 
152  def get_handle_particles(self):
153  return self.handleparticles
154 
155  def get_composite_particles(self):
156  return self.compositeparticles
157 
158  def get_restraint(self):
159  return self.rs
160 
161  def add_to_model(self):
162  IMP.pmi.tools.add_restraint_to_model(self.m, self.rs)
163 
164  def get_output(self):
165  output = {}
166  score = self.rs.unprotected_evaluate(None)
167  output["_TotalScore"] = str(score)
168  output["CompositeRestraint_" + self.label] = str(score)
169  return output
170 
171 
172 #
174 
175  '''
176  this restraint allows ambiguous crosslinking between multiple copies
177  excluding between symmetric copies
178  It allows name ambiguity
179  '''
180 
181  def __init__(
182  self,
183  representation,
184  restraints_file,
185  cut_off=5.0,
186  lam=1.0,
187  plateau=0.01,
188  resolution=None,
189  label="None"):
190 
191  self.weight = 1.0
192  self.m = representation.prot.get_model()
193  self.rs = IMP.RestraintSet(self.m, 'data')
194  self.label = "None"
195  self.pairs = []
196 
197  self.outputlevel = "low"
198  self.cut_off = cut_off
199  self.lam = lam
200  self.plateau = plateau
201 
202  fl = IMP.pmi.tools.open_file_or_inline_text(restraints_file)
203 
204  for line in fl:
205 
206  tokens = line.split()
207  # skip character
208  if (tokens[0] == "#"):
209  continue
210  r1 = int(tokens[2])
211  c1 = tokens[0]
212  r2 = int(tokens[3])
213  c2 = tokens[1]
214 
215  ps1 = IMP.pmi.tools.select(
216  representation,
217  resolution=resolution,
218  name=c1,
219  name_is_ambiguous=True,
220  residue=r1)
221  hrc1 = [representation.hier_db.particle_to_name[p] for p in ps1]
222  ps1nosym = [
223  p for p in ps1 if IMP.pmi.Symmetric(
224  p).get_symmetric(
225  ) == 0]
226  hrc1nosym = [representation.hier_db.particle_to_name[p]
227  for p in ps1nosym]
228 
229  if len(ps1) == 0:
230  print(
231  "AmbiguousCompositeRestraint: WARNING> residue %d of chain %s is not there" %
232  (r1, c1))
233  continue
234 
235  ps2 = IMP.pmi.tools.select(
236  representation,
237  resolution=resolution,
238  name=c2,
239  name_is_ambiguous=True,
240  residue=r2)
241  hrc2 = [representation.hier_db.particle_to_name[p] for p in ps2]
242  ps2nosym = [
243  p for p in ps2 if IMP.pmi.Symmetric(
244  p).get_symmetric(
245  ) == 0]
246  hrc2nosym = [representation.hier_db.particle_to_name[p]
247  for p in ps2nosym]
248 
249  if len(ps2) == 0:
250  print(
251  "AmbiguousCompositeRestraint: WARNING> residue %d of chain %s is not there" %
252  (r2, c2))
253  continue
254 
256  self.m,
257  ps1nosym,
258  self.cut_off,
259  self.lam,
260  True,
261  self.plateau)
262  cr.add_composite_particle(ps2)
263 
264  self.rs.add_restraint(cr)
265  self.pairs.append(
266  (ps1nosym,
267  hrc1nosym,
268  c1,
269  r1,
270  ps2,
271  hrc2,
272  c2,
273  r2,
274  cr))
275 
277  self.m,
278  ps1,
279  self.cut_off,
280  self.lam,
281  True,
282  self.plateau)
283  cr.add_composite_particle(ps2nosym)
284 
285  self.rs.add_restraint(cr)
286  self.pairs.append(
287  (ps1,
288  hrc1,
289  c1,
290  r1,
291  ps2nosym,
292  hrc2nosym,
293  c2,
294  r2,
295  cr))
296 
297  def plot_restraint(
298  self,
299  maxdist=100,
300  npoints=100):
301 
302  p1 = IMP.Particle(self.m)
303  p2 = IMP.Particle(self.m)
307  self.m,
308  [p1],
309  self.cut_off,
310  self.lam,
311  True,
312  self.plateau)
313  cr.add_composite_particle([p2])
314  dists = []
315  scores = []
316  for i in range(npoints):
317  d2.set_coordinates(
318  IMP.algebra.Vector3D(maxdist / npoints * float(i), 0, 0))
319  dists.append(IMP.core.get_distance(d1, d2))
320  scores.append(cr.unprotected_evaluate(None))
321  IMP.pmi.output.plot_xy_data(dists, scores)
322 
323  def set_label(self, label):
324  self.label = label
325  self.rs.set_name(label)
326  for r in self.rs.get_restraints():
327  r.set_name(label)
328 
329  def add_to_model(self):
330  IMP.pmi.tools.add_restraint_to_model(self.m, self.rs)
331 
332  def get_hierarchies(self):
333  return self.prot
334 
335  def get_restraint_sets(self):
336  return self.rs
337 
338  def get_restraint(self):
339  return self.rs
340 
341  def set_output_level(self, level="low"):
342  # this might be "low" or "high"
343  self.outputlevel = level
344 
345  def set_weight(self, weight):
346  self.weight = weight
347  self.rs.set_weight(weight)
348 
349  def get_output(self):
350  # content of the crosslink database pairs
351  # self.pairs.append((p1,p2,dr,r1,c1,r2,c2))
352  output = {}
353  score = self.weight * self.rs.unprotected_evaluate(None)
354  output["_TotalScore"] = str(score)
355  output["AmbiguousCompositeRestraint_Score_" + self.label] = str(score)
356  for n, p in enumerate(self.pairs):
357 
358  ps1 = p[0]
359  hrc1 = p[1]
360  c1 = p[2]
361  r1 = p[3]
362  ps2 = p[4]
363  hrc2 = p[5]
364  c2 = p[6]
365  r2 = p[7]
366  cr = p[8]
367  for n1, p1 in enumerate(ps1):
368  name1 = hrc1[n1]
369 
370  for n2, p2 in enumerate(ps2):
371  name2 = hrc2[n2]
372  d1 = IMP.core.XYZR(p1)
373  d2 = IMP.core.XYZR(p2)
374  label = str(r1) + ":" + name1 + "_" + str(r2) + ":" + name2
375  output["AmbiguousCompositeRestraint_Distance_" +
376  label] = str(IMP.core.get_distance(d1, d2))
377 
378  label = str(r1) + ":" + c1 + "_" + str(r2) + ":" + c2
379  output["AmbiguousCompositeRestraint_Score_" +
380  label] = str(self.weight * cr.unprotected_evaluate(None))
381 
382  return output
383 
384 
385 #
386 class SimplifiedPEMAP(object):
387 
388  def __init__(
389  self,
390  representation,
391  restraints_file,
392  expdistance,
393  strength,
394  resolution=None):
395 
396  self.m = representation.prot.get_model()
397  self.rs = IMP.RestraintSet(self.m, 'data')
398  self.label = "None"
399  self.pairs = []
400 
401  self.outputlevel = "low"
402  self.expdistance = expdistance
403  self.strength = strength
404 
405  fl = IMP.pmi.tools.open_file_or_inline_text(restraints_file)
406 
407  for line in fl:
408 
409  tokens = line.split()
410  # skip character
411  if (tokens[0] == "#"):
412  continue
413  r1 = int(tokens[2])
414  c1 = tokens[0]
415  r2 = int(tokens[3])
416  c2 = tokens[1]
417  pcc = float(tokens[4])
418 
419  ps1 = IMP.pmi.tools.select(
420  representation,
421  resolution=resolution,
422  name=c1,
423  name_is_ambiguous=False,
424  residue=r1)
425  if len(ps1) == 0:
426  print(
427  "SimplifiedPEMAP: WARNING> residue %d of chain %s is not there (w/ %d %s)" %
428  (r1, c1, r2, c2))
429  continue
430  if len(ps1) > 1:
431  print(
432  "SimplifiedPEMAP: WARNING> residue %d of chain %s selected multiple particles" %
433  (r1, c1))
434  continue
435 
436  ps2 = IMP.pmi.tools.select(
437  representation,
438  resolution=resolution,
439  name=c2,
440  name_is_ambiguous=False,
441  residue=r2)
442  if len(ps2) == 0:
443  print(
444  "SimplifiedPEMAP: WARNING> residue %d of chain %s is not there (w/ %d %s)" %
445  (r1, c1, r2, c2))
446  continue
447  if len(ps2) > 1:
448  print(
449  "SimplifiedPEMAP: WARNING> residue %d of chain %s selected multiple particles" %
450  (r2, c2))
451  continue
452 
453  p1 = ps1[0]
454  p2 = ps2[0]
455 
456  # This is harmonic potential for the pE-MAP data
457  upperdist = self.get_upper_bond(pcc)
458  limit = self.strength * (upperdist + 15) ** 2 + 10.0
460  upperdist,
461  self.strength,
462  upperdist +
463  15,
464  limit)
465 
466  # This is harmonic for the X-link
467  #hub= IMP.core.TruncatedHarmonicBound(17.0,self.strength,upperdist+15,limit)
468 
470  dr = IMP.core.PairRestraint(self.m, df, (p1, p2))
471  self.rs.add_restraint(dr)
472  self.pairs.append((p1, p2, dr, r1, c1, r2, c2))
473 
474  # Lower-bound restraint
475  lowerdist = self.get_lower_bond(pcc)
476  limit = self.strength * (lowerdist - 15) ** 2 + 10.0
478  lowerdist,
479  self.strength,
480  lowerdist +
481  15,
482  limit)
483 
484  # This is harmonic for the X-link
485  #hub2= IMP.core.TruncatedHarmonicBound(17.0,self.strength,upperdist+15,limit)
486 
488  dr2 = IMP.core.PairRestraint(df2, (p1, p2))
489  self.rs.add_restraint(dr2)
490  self.pairs.append((p1, p2, dr2, r1, c1, r2, c2))
491 
492  def get_upper_bond(self, pearsoncc):
493  # return (pearsoncc-1.)/-0.0075
494  return (pearsoncc - .5) / (-0.005415)
495 
496  def get_lower_bond(self, pearsoncc):
497  return (pearsoncc - 1.) / -0.0551
498 
499  def set_label(self, label):
500  self.label = label
501 
502  def add_to_model(self):
503  IMP.pmi.tools.add_restraint_to_model(self.m, self.rs)
504 
505  def get_hierarchies(self):
506  return self.prot
507 
508  def get_restraint_sets(self):
509  return self.rs
510 
511  def set_output_level(self, level="low"):
512  # this might be "low" or "high"
513  self.outputlevel = level
514 
515  def get_output(self):
516  # content of the crosslink database pairs
517  # self.pairs.append((p1,p2,dr,r1,c1,r2,c2))
518  output = {}
519  score = self.rs.unprotected_evaluate(None)
520  output["_TotalScore"] = str(score)
521  output["SimplifiedPEMAP_Score_" + self.label] = str(score)
522  for i in range(len(self.pairs)):
523 
524  p0 = self.pairs[i][0]
525  p1 = self.pairs[i][1]
526  crosslinker = 'standard'
527  ln = self.pairs[i][2]
528  resid1 = self.pairs[i][3]
529  chain1 = self.pairs[i][4]
530  resid2 = self.pairs[i][5]
531  chain2 = self.pairs[i][6]
532 
533  label = str(resid1) + ":" + chain1 + "_" + \
534  str(resid2) + ":" + chain2
535  output["SimplifiedPEMAP_Score_" + crosslinker + "_" +
536  label] = str(ln.unprotected_evaluate(None))
537 
538  d0 = IMP.core.XYZ(p0)
539  d1 = IMP.core.XYZ(p1)
540  output["SimplifiedPEMAP_Distance_" +
541  label] = str(IMP.core.get_distance(d0, d1))
542 
543  return output
544 
545 
547 
548  '''
549  generates and wraps a IMP.pmi.ConnectivityRestraint between domains
550  example:
551  cr=restraints.ConnectivityNetworkRestraint(simo,["CCC",(1,100,"TTT"),(100,150,"AAA")])
552  cr.add_to_model()
553  cr.set_label("CR1")
554 
555  Multistate support =No
556  Selection type=selection tuple
557  Resolution=Yes
558  '''
559 
560  def __init__(
561  self,
562  representation=None,
563  objects=None,
564  kappa=10.0,
565  resolution=1.0,
566  label="None"):
567 
568  self.weight = 1.0
569  self.kappa = kappa
570  self.label = label
571  if self.label == "None":
572  self.label = str(selection_tuples)
573 
574  hiers=[]
575 
576  if representation is None:
577  print(objects)
578  for obj in objects:
579  hiers.append(IMP.pmi.tools.input_adaptor(obj,
580  resolution,
581  flatten=True))
582  self.m=hiers[0][0].get_model()
583  else:
584  self.m = representation.m
585  for s in objects:
586  hiers.append(IMP.pmi.tools.select_by_tuple(representation, s,
587  resolution=resolution,
588  name_is_ambiguous=False))
589 
590  #particles=[h.get_particle() for h in hiers]
591  cr = ConnectivityNetworkRestraint(self.m)
592  for hs in hiers:
593  cr.add_particles([h.get_particle() for h in hs])
594  self.rs = IMP.RestraintSet(self.m, label)
595  self.rs.add_restraint(cr)
596 
597  def set_label(self, label):
598  self.label = label
599  self.rs.set_name(label)
600  for r in self.rs.get_restraints():
601  r.set_name(label)
602 
603  def add_to_model(self):
604  IMP.pmi.tools.add_restraint_to_model(self.m, self.rs)
605 
606  def get_restraint(self):
607  return self.rs
608 
609  def get_restraints(self):
610  rlist = []
611  for r in self.rs.get_restraints():
612  rlist.append(IMP.core.PairRestraint.get_from(r))
613  return rlist
614 
615  def set_weight(self, weight):
616  self.weight = weight
617  self.rs.set_weight(weight)
618 
619  def get_output(self):
620  output = {}
621  score = self.weight * self.rs.unprotected_evaluate(None)
622  output["_TotalScore"] = str(score)
623  output["ConnectivityNetworkRestraint_" + self.label] = str(score)
624  return output
625 
626 
628 
629  '''
630  a python restraint that computes the score for a composite of proteins
631  Authors: G. Bouvier, R. Pellarin. Pasteur Institute.
632  '''
633 
634  def __init__(self, m, slope=1.0, theta=0.0, plateau=0.0000000001, linear_slope=0.015):
635  '''
636  input a list of particles, the slope and theta of the sigmoid potential
637  theta is the cutoff distance for a protein-protein contact
638  '''
639  # Import networkx here so that we don't introduce it as a dependency
640  # for *every* proteomics restraint, only this one
641  import networkx
642  self.networkx = networkx
643  IMP.Restraint.__init__(self, m, "ConnectivityNetworkRestraint %1%")
644  self.slope = slope
645  self.theta = theta
646  self.linear_slope = linear_slope
647  self.plateau = plateau
648  self.particles_blocks = []
649  self.particle_list = []
650 
651  def get_number_of_particle_blocks(self):
652  return len(self.particles_blocks)
653 
654  def get_number_of_particles_for_block(self, block_index):
655  return len(self.particles_blocks[block_index])
656 
657  def add_particles(self, particles):
658  self.particles_blocks.append(particles)
659  self.particle_list += particles
660 
661  def get_full_graph(self):
662  '''
663  get the full graph of distances between every particle pair
664  '''
665  import scipy.spatial
666  pdist_array = numpy.array(
668  pdist_mat = scipy.spatial.distance.squareform(pdist_array)
669  pdist_mat[pdist_mat < 0] = 0
670  graph = self.networkx.Graph(pdist_mat)
671  return graph
672 
674  """
675  return the minimum spanning tree
676  """
677  graph = self.get_full_graph()
678  graph = self.networkx.minimum_spanning_tree(graph)
679  return graph
680 
681  def sigmoid(self, x):
682  '''
683  a sigmoid function that scores the probability of a contact
684  between two proteins
685  '''
686  # return 1 - (x)**self.slope/ float(((x)**self.slope +
687  # self.theta**self.slope))
688  argvalue = (x - self.theta) / self.slope
689  return 1.0 - (1.0 - self.plateau) / (1.0 + math.exp(-argvalue))
690 
691  def unprotected_evaluate(self, da):
692  graph = self.get_minimum_spanning_tree()
693  score = 0.0
694  for e in graph.edges():
695  dist = graph.get_edge_data(*e)['weight']
696  prob = self.sigmoid(dist)
697  score += -numpy.log(prob)
698  score += self.linear_slope * dist
699  return score
700 
701  def do_get_inputs(self):
702  return self.particle_list
703 
705  "use IMP.restraints.basic.MembraneRestraint instead")
706 class SetupMembraneRestraint(object):
707 
708  '''
709 
710  '''
711 
712  def get_from_selection_tuple(self,tuples):
713  particles = []
714  for s in tuples:
715  ps = IMP.pmi.tools.select_by_tuple(
716  self.representation, s,
717  resolution=self.resolution, name_is_ambiguous=True)
718  particles += ps
719  return particles
720 
721  def __init__(
722  self,
723  representation=None,
724  objects_above=None,
725  objects_inside=None,
726  objects_below=None,
727  z_init=0.0,
728  z_min=0.0,
729  z_max=0.0,
730  thickness=30,
731  resolution=1,
732  label="None"):
733  self.resolution=resolution
734 
735  self.weight = 1.0
736  self.label = label
737  self.representation = representation
738  self.thickness = thickness
739  self.z_init=z_init
740 
741 
742  softness = 3.0
743  plateau = 1e-10
744  linear = 0.02
745 
746  mr=None
747  hierarchies_above = []
748  hierarchies_inside = []
749  hierarchies_below = []
750 
751  if representation is None:
752  hierarchies_above, hierarchies_inside, hierarchies_below = (
753  IMP.pmi.tools.input_adaptor(objects,
754  resolution,
755  flatten=True) for objects in [objects_above, objects_inside, objects_below])
756  for h in hierarchies_above, hierarchies_inside, hierarchies_below:
757  if h is not None:
758  self.m=hierarchies_above[0].get_model()
759  break
760  else:
761  self.m = representation.prot.get_model()
762  if objects_above is not None:
763  hierarchies_above = self.get_from_selection_tuple(objects_above)
764 
765  if objects_below is not None:
766  particles_below = self.get_from_selection_tuple(objects_below)
767 
768  if objects_inside is not None:
769  particles_inside = self.get_from_selection_tuple(objects_inside)
770 
771  self.z_center = IMP.pmi.tools.SetupNuisance(
772  self.m, self.z_init, z_min, z_max, isoptimized=True).get_particle()
774  self.m, self.z_center.get_particle_index(), self.thickness, softness, plateau, linear)
775  mr.add_particles_inside([h.get_particle().get_index()
776  for h in hierarchies_inside])
777  mr.add_particles_above([h.get_particle().get_index()
778  for h in hierarchies_above])
779  mr.add_particles_below([h.get_particle().get_index()
780  for h in hierarchies_below])
781 
782  self.rs = IMP.RestraintSet(self.m, label)
783  self.rs.add_restraint(mr)
784 
785  def create_box(self, x_center, y_center, hierarchy):
786 
787  z = self.z_center.get_nuisance()
788  p = IMP.Particle(self.m)
790  h.set_name("Membrane_" + self.label)
791  hierarchy.add_child(h)
792 
793  particles_box = []
794  p_origin = IMP.Particle(self.m)
795 
796  IMP.atom.Mass.setup_particle(p_origin, 100)
797  d = IMP.core.XYZR.setup_particle(p_origin)
798  d.set_coordinates((x_center, y_center, 0))
799  d.set_radius(1)
800  h_origin = IMP.atom.Hierarchy.setup_particle(p_origin)
801  h.add_child(h_origin)
802  particles_box.append(p_origin)
803 
804  p1 = IMP.Particle(self.m)
807  d.set_coordinates((x_center, y_center, z + 100))
808  d.set_radius(1)
810  h.add_child(h1)
811  particles_box.append(p1)
812 
813  p2 = IMP.Particle(self.m)
816  d.set_coordinates((x_center, y_center, z - 100))
817  d.set_radius(1)
819  h.add_child(h2)
820  particles_box.append(p2)
821 
824 
826 
827  p1.set_name("z_positive")
828  p2.set_name("z_negative")
829 
830  # to display the membrane in the rmf file
831  for offs in [self.thickness / 2, -self.thickness / 2]:
832  p1 = IMP.Particle(self.m)
835  d.set_coordinates((-100 + x_center, -100 + y_center, z + offs))
836  d.set_radius(1)
838  h.add_child(h1)
839  particles_box.append(p1)
840 
841  p2 = IMP.Particle(self.m)
844  d.set_coordinates((-100 + x_center, 100 + y_center, z + offs))
845  d.set_radius(1)
847  h.add_child(h2)
848  particles_box.append(p2)
849 
850  p3 = IMP.Particle(self.m)
853  d.set_coordinates((100 + x_center, -100 + y_center, z + offs))
854  d.set_radius(1)
856  h.add_child(h3)
857  particles_box.append(p3)
858 
859  p4 = IMP.Particle(self.m)
862  d.set_coordinates((100 + x_center, 100 + y_center, z + offs))
863  d.set_radius(1)
865  h.add_child(h4)
866  particles_box.append(p4)
867 
872 
877 
878  #sm = self._MembraneSingletonModifier(p_origin, self.z_center)
880  for p in particles_box:
881  smp = self._MembraneSingletonModifier(IMP.core.XYZ(p).get_z(),self.z_init, self.z_center)
883  IMP.core.XYZ(p).set_coordinates_are_optimized(True)
884  lcp.add(p.get_index())
885  c = IMP.container.SingletonsConstraint(smp, None, lcp)
886  self.m.add_score_state(c)
887  self.m.update()
888 
889  class _MembraneSingletonModifier(IMP.SingletonModifier):
890 
891  """A class that updates the membrane particles
892  """
893 
894  def __init__(self,z, z_init, z_nuisance):
895  IMP.SingletonModifier.__init__(
896  self, "MembraneSingletonModifier%1%")
897  self.z=z
898  self.z_init = z_init
899  self.z_nuisance = z_nuisance
900 
901  def apply_index(self, m, pi):
902  new_z = self.z_nuisance.get_nuisance()
903  d = IMP.core.XYZ(m, pi)
904  d.set_z(self.z + new_z-self.z_init)
905 
906 
907 
908  def do_get_inputs(self, m, pis):
909  return IMP.get_particles(m, pis)
910 
911  def do_get_outputs(self, m, pis):
912  return self.do_get_inputs(m, pis)
913 
914  def set_label(self, label):
915  self.label = label
916  self.rs.set_name(label)
917  for r in self.rs.get_restraints():
918  r.set_name(label)
919 
920  def add_to_model(self):
921  IMP.pmi.tools.add_restraint_to_model(self.m, self.rs)
922 
923  def get_restraint(self):
924  return self.rs
925 
926  def set_weight(self, weight):
927  self.weight = weight
928  self.rs.set_weight(weight)
929 
930  def get_particles_to_sample(self):
931  ps = {}
932 
933  ps["Nuisances_MembraneRestraint_Z_" +
934  self.label] = ([self.z_center], 2.0)
935  return ps
936 
937  def get_movers(self):
938  movers=[]
939  mover_name="Nuisances_MembraneRestraint_Z_" + self.label
940  particle=self.z_center
941  maxstep=2.0
942  mv=IMP.core.NormalMover([particle],
943  IMP.FloatKeys([IMP.FloatKey("nuisance")]),maxstep)
944  mv.set_name(mover_name)
945  movers.append(mv)
946 
947  return movers
948 
949 
950  def get_output(self):
951  output = {}
952  score = self.weight * self.rs.unprotected_evaluate(None)
953  output["_TotalScore"] = str(score)
954  output["MembraneRestraint_" + self.label] = str(score)
955  output["MembraneRestraint_Z_" +
956  self.label] = str(self.z_center.get_nuisance())
957  return output
958 
959 
960 class FuzzyBoolean(object):
961 
962  '''
963  Fully Ambiguous Restraint that can be built using boolean logic
964  R. Pellarin. Pasteur Institute.
965  '''
966 
967  def __init__(self, p1, operator=None, p2=None):
968  '''
969  input a list of particles, the slope and theta of the sigmoid potential
970  theta is the cutoff distance for a protein-protein contact
971  '''
972  if isinstance(p1, FuzzyBoolean) and isinstance(p2, FuzzyBoolean):
973  self.operations = [p1, operator, p2]
974  self.value = None
975  else:
976  self.operations = []
977  self.value = p1
978 
979  def __or__(self, FuzzyBoolean2):
980  return FuzzyBoolean(self, self.or_, FuzzyBoolean2)
981 
982  def __and__(self, FuzzyBoolean2):
983  return FuzzyBoolean(self, self.and_, FuzzyBoolean2)
984 
985  def and_(self, a, b):
986  return a * b
987 
988  def or_(self, a, b):
989  return 1.0 - (1.0 - a) * (1.0 - b)
990 
991  def evaluate(self):
992 
993  if len(self.operations) == 0:
994  return self.value
995  FuzzyBoolean1, op, FuzzyBoolean2 = self.operations
996 
997  return op(FuzzyBoolean1.evaluate(), FuzzyBoolean2.evaluate())
998 
999 
1001 
1002  '''
1003  Fully Ambiguous Restraint that can be built using boolean logic
1004  R. Pellarin. Pasteur Institute.
1005  '''
1006  plateau = 0.00000000000001
1007  theta = 5.0
1008  slope = 2.0
1009  innerslope = 0.01
1010 
1011  def __init__(self, m, p1, p2, operator=None):
1012  '''
1013  input a list of particles, the slope and theta of the sigmoid potential
1014  theta is the cutoff distance for a protein-protein contact
1015  '''
1016  IMP.Restraint.__init__(self, m, "FuzzyRestraint %1%")
1017  self.m = m
1018  self.min = sys.float_info.min
1019  if isinstance(p1, FuzzyRestraint) and isinstance(p2, FuzzyRestraint):
1020  self.operations = [p1, operator, p2]
1021  self.particle_pair = None
1022  elif isinstance(p1, FuzzyRestraint) and p2 is None:
1023  self.operations = [p1, operator, None]
1024  self.particle_pair = None
1025  else:
1026  self.operations = []
1027  self.particle_pair = (p1, p2)
1028 
1029  def __or__(self, FuzzyRestraint2):
1030  return FuzzyRestraint(self.m, self, FuzzyRestraint2, self.or_)
1031 
1032  def __and__(self, FuzzyRestraint2):
1033  return FuzzyRestraint(self.m, self, FuzzyRestraint2, self.and_)
1034 
1035  def __invert__(self):
1036  return FuzzyRestraint(self.m, self, None, self.invert_)
1037 
1038  def and_(self, a, b):
1039  c = a + b
1040  return c
1041 
1042  def or_(self, a, b):
1043  c = math.exp(-a) + math.exp(-b) - math.exp(-a - b)
1044  return -math.log(c)
1045 
1046  def invert_(self, a):
1047  c = 1.0 - math.exp(-a)
1048  return -math.log(c)
1049 
1050  def evaluate(self):
1051  if len(self.operations) == 0:
1052  return self.distance()
1053  FuzzyRestraint1, op, FuzzyRestraint2 = self.operations
1054 
1055  if FuzzyRestraint2 is not None:
1056  return op(FuzzyRestraint1.evaluate(), FuzzyRestraint2.evaluate())
1057  else:
1058  return op(FuzzyRestraint1.evaluate())
1059 
1060  def distance(self):
1061  d1 = IMP.core.XYZ(self.particle_pair[0])
1062  d2 = IMP.core.XYZ(self.particle_pair[1])
1063  d = IMP.core.get_distance(d1, d2)
1064  argvalue = (d-self.theta)/self.slope
1065  return -math.log(1.0 -(1.0-self.plateau)/(1.0+math.exp(-argvalue)))+self.innerslope*d
1066 
1067  def add_to_model(self):
1069 
1070  def unprotected_evaluate(self, da):
1071  return self.evaluate()
1072 
1073  def __str__(self):
1074  if len(self.operations) == 0:
1075  return str(self.particle_pair)
1076  FuzzyRestraint1, op, FuzzyRestraint2 = self.operations
1077  if FuzzyRestraint2 is not None:
1078  return str(FuzzyRestraint1) +str(op)+str(FuzzyRestraint2)
1079  else:
1080  return str(FuzzyRestraint1) +str(op)
1081 
1082  def do_get_inputs(self):
1083  if len(self.operations) == 0:
1084  return list(self.particle_pair)
1085  FuzzyRestraint1, op, FuzzyRestraint2 = self.operations
1086  if FuzzyRestraint2 is not None:
1087  return list(set(FuzzyRestraint1.do_get_inputs() +FuzzyRestraint2.do_get_inputs()))
1088  else:
1089  return list(set(FuzzyRestraint1.do_get_inputs()))
A base class for modifiers of ParticlesTemp.
A function that is harmonic over an interval.
Apply a SingletonFunction to a SingletonContainer to maintain an invariant.
this restraint allows ambiguous crosslinking between multiple copies excluding between symmetric copi...
void add_particles(RMF::FileHandle fh, const ParticlesTemp &hs)
Various classes to hold sets of particles.
static XYZR setup_particle(Model *m, ParticleIndex pi)
Definition: XYZR.h:48
A decorator for a particle which has bonds.
generate a connectivity restraint between domains setting up the restraint example: cr=restraints...
def get_full_graph
get the full graph of distances between every particle pair
Miscellaneous utilities.
Definition: tools.py:1
Add symmetric attribute to a particle.
Definition: Symmetric.h:23
ParticlesTemp get_particles(Model *m, const ParticleIndexes &ps)
A score on the distance between the surfaces of two spheres.
double get_distance(XYZR a, XYZR b)
Compute the sphere distance between a and b.
Definition: XYZR.h:89
generates and wraps a IMP.pmi.ConnectivityRestraint between domains example: cr=restraints.ConnectivityNetworkRestraint(simo,["CCC",(1,100,"TTT"),(100,150,"AAA")]) cr.add_to_model() cr.set_label("CR1")
Bond create_bond(Bonded a, Bonded b, Bond o)
Connect the two wrapped particles by a custom bond.
Object used to hold a set of restraints.
Definition: RestraintSet.h:36
a python restraint that computes the score for a composite of proteins Authors: G.
static Hierarchy setup_particle(Model *m, ParticleIndex pi, ParticleIndexesAdaptor children=ParticleIndexesAdaptor())
Create a Hierarchy of level t by adding the needed attributes.
def deprecated_pmi1_object
Mark a PMI1 class as deprecated.
Ints get_index(const ParticlesTemp &particles, const Subset &subset, const Subsets &excluded)
def __init__
input a list of particles, the slope and theta of the sigmoid potential theta is the cutoff distance ...
Store a list of ParticleIndexes.
static Mass setup_particle(Model *m, ParticleIndex pi, Float mass)
Definition: Mass.h:44
def add_restraint_to_model
Add a PMI restraint to the model.
Definition: tools.py:69
static Bonded setup_particle(Model *m, ParticleIndex pi)
A restraint for ambiguous cross-linking MS data and multiple state approach.
A decorator for a particle with x,y,z coordinates.
Definition: XYZ.h:30
handleparticles is a selection tuple compositeparticles is a list of selection tuples ...
Modify a set of continuous variables using a normal distribution.
Definition: NormalMover.h:20
Classes for writing output files and processing them.
Definition: output.py:1
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...
Fully Ambiguous Restraint that can be built using boolean logic R.
def __init__
input a list of particles, the slope and theta of the sigmoid potential theta is the cutoff distance ...
VectorD< 3 > Vector3D
Definition: VectorD.h:395
Class to handle individual particles of a Model object.
Definition: Particle.h:41
Floats get_list_of_bipartite_minimum_sphere_distance(const ParticlesTemps &pss)
Definition: utilities.h:71
def sigmoid
a sigmoid function that scores the probability of a contact between two proteins
Restraint * create_connectivity_restraint(const Selections &s, double x0, double k, std::string name="Connectivity%1%")
Create a restraint connecting the selections.
Applies a PairScore to a Pair.
Definition: PairRestraint.h:29
def select
this function uses representation=SimplifiedModel it returns the corresponding selected particles rep...
Definition: tools.py:726
Python classes to represent, score, sample and analyze models.
Functionality for loading, creating, manipulating and scoring atomic structures.
Select hierarchy particles identified by the biological name.
Definition: Selection.h:66
Fully Ambiguous Restraint that can be built using boolean logic R.
def __init__
input a list of particles, the slope and theta of the sigmoid potential theta is the cutoff distance ...
virtual ModelObjectsTemp do_get_inputs() const =0
A restraint is a term in an IMP ScoringFunction.
Definition: Restraint.h:54
A decorator for a particle with x,y,z coordinates and a radius.
Definition: XYZR.h:27