IMP logo
IMP Reference Guide  2.10.0
The Integrative Modeling Platform
pmi1/restraints/proteomics.py
1 """@namespace IMP.pmi1.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.pmi1
12 import IMP.pmi1.tools
13 import IMP.pmi1.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.pmi1.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.pmi1.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.pmi1.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):
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.pmi1.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.pmi1.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.pmi1.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.pmi1.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.pmi1.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.pmi1.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):
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.pmi1.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.pmi1.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.pmi1.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(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):
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.pmi1.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.pmi1.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.pmi1.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):
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 
704 
706 
707  '''
708 
709  '''
710 
711  def get_from_selection_tuple(self,tuples):
712  particles = []
713  for s in tuples:
714  ps = IMP.pmi1.tools.select_by_tuple(
715  self.representation, s,
716  resolution=self.resolution, name_is_ambiguous=True)
717  particles += ps
718  return particles
719 
720  def __init__(
721  self,
722  representation=None,
723  objects_above=None,
724  objects_inside=None,
725  objects_below=None,
726  z_init=0.0,
727  z_min=0.0,
728  z_max=0.0,
729  thickness=30,
730  resolution=1,
731  label="None"):
732  self.resolution=resolution
733 
734  self.weight = 1.0
735  self.label = label
736  self.representation = representation
737  self.thickness = thickness
738  self.z_init=z_init
739 
740 
741  softness = 3.0
742  plateau = 1e-10
743  linear = 0.02
744 
745  mr=None
746  hierarchies_above = []
747  hierarchies_inside = []
748  hierarchies_below = []
749 
750  if representation is None:
751  hierarchies_above, hierarchies_inside, hierarchies_below = (
753  resolution,
754  flatten=True) for objects in [objects_above, objects_inside, objects_below])
755  for h in hierarchies_above, hierarchies_inside, hierarchies_below:
756  if h is not None:
757  self.m=hierarchies_above[0].get_model()
758  break
759  else:
760  self.m = representation.prot.get_model()
761  if objects_above is not None:
762  hierarchies_above = self.get_from_selection_tuple(objects_above)
763 
764  if objects_below is not None:
765  particles_below = self.get_from_selection_tuple(objects_below)
766 
767  if objects_inside is not None:
768  particles_inside = self.get_from_selection_tuple(objects_inside)
769 
770  self.z_center = IMP.pmi1.tools.SetupNuisance(
771  self.m, self.z_init, z_min, z_max, isoptimized=True).get_particle()
773  self.m, self.z_center.get_particle_index(), self.thickness, softness, plateau, linear)
774  mr.add_particles_inside([h.get_particle().get_index()
775  for h in hierarchies_inside])
776  mr.add_particles_above([h.get_particle().get_index()
777  for h in hierarchies_above])
778  mr.add_particles_below([h.get_particle().get_index()
779  for h in hierarchies_below])
780 
781  self.rs = IMP.RestraintSet(self.m, label)
782  self.rs.add_restraint(mr)
783 
784  def create_box(self, x_center, y_center, hierarchy):
785 
786  z = self.z_center.get_nuisance()
787  p = IMP.Particle(self.m)
789  h.set_name("Membrane_" + self.label)
790  hierarchy.add_child(h)
791 
792  particles_box = []
793  p_origin = IMP.Particle(self.m)
794 
795  IMP.atom.Mass.setup_particle(p_origin, 100)
796  d = IMP.core.XYZR.setup_particle(p_origin)
797  d.set_coordinates((x_center, y_center, 0))
798  d.set_radius(1)
799  h_origin = IMP.atom.Hierarchy.setup_particle(p_origin)
800  h.add_child(h_origin)
801  particles_box.append(p_origin)
802 
803  p1 = IMP.Particle(self.m)
806  d.set_coordinates((x_center, y_center, z + 100))
807  d.set_radius(1)
809  h.add_child(h1)
810  particles_box.append(p1)
811 
812  p2 = IMP.Particle(self.m)
815  d.set_coordinates((x_center, y_center, z - 100))
816  d.set_radius(1)
818  h.add_child(h2)
819  particles_box.append(p2)
820 
823 
825 
826  p1.set_name("z_positive")
827  p2.set_name("z_negative")
828 
829  # to display the membrane in the rmf file
830  for offs in [self.thickness / 2, -self.thickness / 2]:
831  p1 = IMP.Particle(self.m)
834  d.set_coordinates((-100 + x_center, -100 + y_center, z + offs))
835  d.set_radius(1)
837  h.add_child(h1)
838  particles_box.append(p1)
839 
840  p2 = IMP.Particle(self.m)
843  d.set_coordinates((-100 + x_center, 100 + y_center, z + offs))
844  d.set_radius(1)
846  h.add_child(h2)
847  particles_box.append(p2)
848 
849  p3 = IMP.Particle(self.m)
852  d.set_coordinates((100 + x_center, -100 + y_center, z + offs))
853  d.set_radius(1)
855  h.add_child(h3)
856  particles_box.append(p3)
857 
858  p4 = IMP.Particle(self.m)
861  d.set_coordinates((100 + x_center, 100 + y_center, z + offs))
862  d.set_radius(1)
864  h.add_child(h4)
865  particles_box.append(p4)
866 
871 
876 
877  #sm = self._MembraneSingletonModifier(p_origin, self.z_center)
879  for p in particles_box:
880  smp = self._MembraneSingletonModifier(IMP.core.XYZ(p).get_z(),self.z_init, self.z_center)
882  IMP.core.XYZ(p).set_coordinates_are_optimized(True)
883  lcp.add(p.get_index())
884  c = IMP.container.SingletonsConstraint(smp, None, lcp)
885  self.m.add_score_state(c)
886  self.m.update()
887 
888  class _MembraneSingletonModifier(IMP.SingletonModifier):
889 
890  """A class that updates the membrane particles
891  """
892 
893  def __init__(self,z, z_init, z_nuisance):
894  IMP.SingletonModifier.__init__(
895  self, "MembraneSingletonModifier%1%")
896  self.z=z
897  self.z_init = z_init
898  self.z_nuisance = z_nuisance
899 
900  def apply_index(self, m, pi):
901  new_z = self.z_nuisance.get_nuisance()
902  d = IMP.core.XYZ(m, pi)
903  d.set_z(self.z + new_z-self.z_init)
904 
905 
906 
907  def do_get_inputs(self, m, pis):
908  return IMP.get_particles(m, pis)
909 
910  def do_get_outputs(self, m, pis):
911  return self.do_get_inputs(m, pis)
912 
913  def set_label(self, label):
914  self.label = label
915  self.rs.set_name(label)
916  for r in self.rs.get_restraints():
917  r.set_name(label)
918 
919  def add_to_model(self):
921 
922  def get_restraint(self):
923  return self.rs
924 
925  def set_weight(self, weight):
926  self.weight = weight
927  self.rs.set_weight(weight)
928 
929  def get_particles_to_sample(self):
930  ps = {}
931 
932  ps["Nuisances_MembraneRestraint_Z_" +
933  self.label] = ([self.z_center], 2.0)
934  return ps
935 
936  def get_movers(self):
937  movers=[]
938  mover_name="Nuisances_MembraneRestraint_Z_" + self.label
939  particle=self.z_center
940  maxstep=2.0
941  mv=IMP.core.NormalMover([particle],
942  IMP.FloatKeys([IMP.FloatKey("nuisance")]),maxstep)
943  mv.set_name(mover_name)
944  movers.append(mv)
945 
946  return movers
947 
948 
949  def get_output(self):
950  output = {}
951  score = self.weight * self.rs.unprotected_evaluate(None)
952  output["_TotalScore"] = str(score)
953  output["MembraneRestraint_" + self.label] = str(score)
954  output["MembraneRestraint_Z_" +
955  self.label] = str(self.z_center.get_nuisance())
956  return output
957 
958 
959 class FuzzyBoolean(object):
960 
961  '''
962  Fully Ambiguous Restraint that can be built using boolean logic
963  R. Pellarin. Pasteur Institute.
964  '''
965 
966  def __init__(self, p1, operator=None, p2=None):
967  '''
968  input a list of particles, the slope and theta of the sigmoid potential
969  theta is the cutoff distance for a protein-protein contact
970  '''
971  if isinstance(p1, FuzzyBoolean) and isinstance(p2, FuzzyBoolean):
972  self.operations = [p1, operator, p2]
973  self.value = None
974  else:
975  self.operations = []
976  self.value = p1
977 
978  def __or__(self, FuzzyBoolean2):
979  return FuzzyBoolean(self, self.or_, FuzzyBoolean2)
980 
981  def __and__(self, FuzzyBoolean2):
982  return FuzzyBoolean(self, self.and_, FuzzyBoolean2)
983 
984  def and_(self, a, b):
985  return a * b
986 
987  def or_(self, a, b):
988  return 1.0 - (1.0 - a) * (1.0 - b)
989 
990  def evaluate(self):
991 
992  if len(self.operations) == 0:
993  return self.value
994  FuzzyBoolean1, op, FuzzyBoolean2 = self.operations
995 
996  return op(FuzzyBoolean1.evaluate(), FuzzyBoolean2.evaluate())
997 
998 
1000 
1001  '''
1002  Fully Ambiguous Restraint that can be built using boolean logic
1003  R. Pellarin. Pasteur Institute.
1004  '''
1005  plateau = 0.00000000000001
1006  theta = 5.0
1007  slope = 2.0
1008  innerslope = 0.01
1009 
1010  def __init__(self, m, p1, p2, operator=None):
1011  '''
1012  input a list of particles, the slope and theta of the sigmoid potential
1013  theta is the cutoff distance for a protein-protein contact
1014  '''
1015  IMP.Restraint.__init__(self, m, "FuzzyRestraint %1%")
1016  self.m = m
1017  self.min = sys.float_info.min
1018  if isinstance(p1, FuzzyRestraint) and isinstance(p2, FuzzyRestraint):
1019  self.operations = [p1, operator, p2]
1020  self.particle_pair = None
1021  elif isinstance(p1, FuzzyRestraint) and p2 is None:
1022  self.operations = [p1, operator, None]
1023  self.particle_pair = None
1024  else:
1025  self.operations = []
1026  self.particle_pair = (p1, p2)
1027 
1028  def __or__(self, FuzzyRestraint2):
1029  return FuzzyRestraint(self.m, self, FuzzyRestraint2, self.or_)
1030 
1031  def __and__(self, FuzzyRestraint2):
1032  return FuzzyRestraint(self.m, self, FuzzyRestraint2, self.and_)
1033 
1034  def __invert__(self):
1035  return FuzzyRestraint(self.m, self, None, self.invert_)
1036 
1037  def and_(self, a, b):
1038  c = a + b
1039  return c
1040 
1041  def or_(self, a, b):
1042  c = math.exp(-a) + math.exp(-b) - math.exp(-a - b)
1043  return -math.log(c)
1044 
1045  def invert_(self, a):
1046  c = 1.0 - math.exp(-a)
1047  return -math.log(c)
1048 
1049  def evaluate(self):
1050  if len(self.operations) == 0:
1051  return self.distance()
1052  FuzzyRestraint1, op, FuzzyRestraint2 = self.operations
1053 
1054  if FuzzyRestraint2 is not None:
1055  return op(FuzzyRestraint1.evaluate(), FuzzyRestraint2.evaluate())
1056  else:
1057  return op(FuzzyRestraint1.evaluate())
1058 
1059  def distance(self):
1060  d1 = IMP.core.XYZ(self.particle_pair[0])
1061  d2 = IMP.core.XYZ(self.particle_pair[1])
1062  d = IMP.core.get_distance(d1, d2)
1063  argvalue = (d-self.theta)/self.slope
1064  return -math.log(1.0 -(1.0-self.plateau)/(1.0+math.exp(-argvalue)))+self.innerslope*d
1065 
1066  def add_to_model(self):
1068 
1069  def unprotected_evaluate(self, da):
1070  return self.evaluate()
1071 
1072  def __str__(self):
1073  if len(self.operations) == 0:
1074  return str(self.particle_pair)
1075  FuzzyRestraint1, op, FuzzyRestraint2 = self.operations
1076  if FuzzyRestraint2 is not None:
1077  return str(FuzzyRestraint1) +str(op)+str(FuzzyRestraint2)
1078  else:
1079  return str(FuzzyRestraint1) +str(op)
1080 
1081  def do_get_inputs(self):
1082  if len(self.operations) == 0:
1083  return list(self.particle_pair)
1084  FuzzyRestraint1, op, FuzzyRestraint2 = self.operations
1085  if FuzzyRestraint2 is not None:
1086  return list(set(FuzzyRestraint1.do_get_inputs() +FuzzyRestraint2.do_get_inputs()))
1087  else:
1088  return list(set(FuzzyRestraint1.do_get_inputs()))
A base class for modifiers of ParticlesTemp.
Fully Ambiguous Restraint that can be built using boolean logic R.
A function that is harmonic over an interval.
generates and wraps a IMP.pmi1.ConnectivityRestraint between domains example: cr=restraints.ConnectivityNetworkRestraint(simo,["CCC",(1,100,"TTT"),(100,150,"AAA")]) cr.add_to_model() cr.set_label("CR1")
def input_adaptor
Adapt things for PMI (degrees of freedom, restraints, ...) Returns list of list of hierarchies...
Definition: /tools.py:1618
def __init__
input a list of particles, the slope and theta of the sigmoid potential theta is the cutoff distance ...
Apply a SingletonFunction to a SingletonContainer to maintain an invariant.
Floats get_list_of_bipartite_minimum_sphere_distance(const ParticlesTemps &pss)
Definition: /utilities.h:71
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.
Legacy PMI1 module to represent, score, sample and analyze models.
def get_full_graph
get the full graph of distances between every particle pair
ParticlesTemp get_particles(Model *m, const ParticleIndexes &ps)
A score on the distance between the surfaces of two spheres.
def sigmoid
a sigmoid function that scores the probability of a contact between two proteins
double get_distance(XYZR a, XYZR b)
Compute the sphere distance between a and b.
Definition: XYZR.h:89
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.
Miscellaneous utilities.
Definition: /tools.py:1
static Hierarchy setup_particle(Model *m, ParticleIndex pi, ParticleIndexesAdaptor children=ParticleIndexesAdaptor())
Create a Hierarchy of level t by adding the needed attributes.
Ints get_index(const ParticlesTemp &particles, const Subset &subset, const Subsets &excluded)
Store a list of ParticleIndexes.
generate a connectivity restraint between domains setting up the restraint example: cr=restraints...
static Mass setup_particle(Model *m, ParticleIndex pi, Float mass)
Definition: Mass.h:44
Fully Ambiguous Restraint that can be built using boolean logic R.
static Bonded setup_particle(Model *m, ParticleIndex pi)
handleparticles is a selection tuple compositeparticles is a list of selection tuples ...
A decorator for a particle with x,y,z coordinates.
Definition: XYZ.h:30
def __init__
input a list of particles, the slope and theta of the sigmoid potential theta is the cutoff distance ...
A restraint for ambiguous cross-linking MS data and multiple state approach.
Modify a set of continuous variables using a normal distribution.
Definition: NormalMover.h:20
this restraint allows ambiguous crosslinking between multiple copies excluding between symmetric copi...
def select
this function uses representation=SimplifiedModel it returns the corresponding selected particles rep...
Definition: /tools.py:724
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...
VectorD< 3 > Vector3D
Definition: VectorD.h:395
Class to handle individual particles of a Model object.
Definition: Particle.h:41
Add symmetric attribute to a particle.
Definition: /Symmetric.h:23
def __init__
input a list of particles, the slope and theta of the sigmoid potential theta is the cutoff distance ...
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
Classes for writing output files and processing them.
Definition: /output.py:1
def add_restraint_to_model
Add a PMI restraint to the model.
Definition: /tools.py:50
Functionality for loading, creating, manipulating and scoring atomic structures.
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
A decorator for a particle with x,y,z coordinates and a radius.
Definition: XYZR.h:27