IMP logo
IMP Reference Guide  2.21.0
The Integrative Modeling Platform
version 2.21.0
display/displaying_ensembles.py

The script shows a couple experiments with trying to visualize an ensemble of structures. The ensemble is fairly tight on the assembly scale, but there is significant variation between the location and orientation of the individual proteins (which were modeled as rigid bodies). To save space, the models have had their sidechain atoms removed.

1 ## \example display/displaying_ensembles.py
2 # The script shows a couple experiments with trying to visualize an
3 # ensemble of structures. The ensemble is fairly tight on the assembly
4 # scale, but there is significant variation between the location and
5 # orientation of the individual proteins (which were modeled as rigid
6 # bodies). To save space, the models have had their sidechain atoms
7 # removed.
8 
9 from __future__ import print_function
10 import IMP.display
11 import IMP.atom
12 import sys
13 
14 IMP.setup_from_argv(
15  sys.argv, "Experiments with trying to visualize an ensemble of structures")
16 
17 Segment = IMP.algebra.Segment3D
18 Cylinder = IMP.algebra.Cylinder3D
19 
20 # turn off internal checks to speed things up
21 IMP.set_check_level(IMP.USAGE)
22 
23 
24 def read(m, beyond_file):
25  print("reading")
26  hs = []
27  for i in range(0, beyond_file):
28  # create a simplified version for each chain to speed up computations
29  name = IMP.display.get_example_path(
30  "ensemble/aligned-" + str(i) + ".pdb")
31  h = IMP.atom.read_pdb(name, m, IMP.atom.CAlphaPDBSelector())
32  hr = IMP.atom.Hierarchy.setup_particle(IMP.Particle(m))
33  hs.append(hr)
34  for c in IMP.atom.get_by_type(h, IMP.atom.CHAIN_TYPE):
35  simp = IMP.atom.create_simplified_along_backbone(
36  IMP.atom.Chain(c), 4)
37  hr.add_child(simp)
38  IMP.atom.destroy(h)
39  print(" ", i)
40  return hs
41 
42 
43 def add_markers(h, c, w):
44  """Add markers to a the passed conformation. The marker locations are chosen
45  pretty thoughtlessly and don't really illustrate the technique well."""
46  def add_marker(s, name):
47  g = IMP.core.XYZRGeometry(s.get_selected_particles()[0])
48  g.set_name(name)
49  g.set_color(c)
50  w.add_geometry(g)
51  s = IMP.atom.Selection(h, chain='B', residue_index=317)
52  add_marker(s, "m0")
53  s = IMP.atom.Selection(h, chain='G', residue_index=212)
54  add_marker(s, "m1")
55  s = IMP.atom.Selection(h, chain='I', residue_index=237)
56  add_marker(s, "m2")
57  s = IMP.atom.Selection(h, chain='F', residue_index=101)
58  add_marker(s, "m3")
59 
60 
61 def get_nice_name(h):
62  nm = h.get_name()
63  return nm[nm.find('-') + 1:nm.rfind('.')]
64 
65 
66 def add_axis(h, c, w, chain_colors):
67  """Add a coordinate axis to show the relative orientation of the protein"""
68  for hc in IMP.atom.get_by_type(h, IMP.atom.CHAIN_TYPE):
69  rb = IMP.core.RigidMember(hc).get_rigid_body()
70  g = IMP.display.ReferenceFrameGeometry(rb.get_reference_frame())
71  g.set_name(get_nice_name(h) + "_orient")
72  if c:
73  g.set_color(c)
74  else:
75  g.set_color(chain_colors[IMP.atom.Chain(hc).get_id()])
76  w.add_geometry(g)
77 
78 
79 def add_skeleton(h, c, r, w, chain_colors):
80  """Show the connectivity skeleton of the conformation to give an idea of
81  how things are laid out"""
82  for hc0 in IMP.atom.get_by_type(h, IMP.atom.CHAIN_TYPE):
83  for hc1 in IMP.atom.get_by_type(h, IMP.atom.CHAIN_TYPE):
84  if hc1 <= hc0:
85  continue
86  d = ps.evaluate_index(h.get_model(),
87  (hc0.get_particle_index(),
88  hc1.get_particle_index()), None)
89  if d < 1:
90  d0 = IMP.core.XYZ(hc0)
91  d1 = IMP.core.XYZ(hc1)
92  mp = .5 * (d0.get_coordinates() + d1.get_coordinates())
93  g = IMP.display.CylinderGeometry(
94  Cylinder(Segment(d0.get_coordinates(), mp), r))
95  if c:
96  g.set_color(c)
97  else:
98  g.set_color(chain_colors[IMP.atom.Chain(d0).get_id()])
99  g.set_name(get_nice_name(h) + "_skel")
100  w.add_geometry(g)
101  g = IMP.display.CylinderGeometry(
102  Cylinder(Segment(d1.get_coordinates(), mp), r))
103  if c:
104  g.set_color(c)
105  else:
106  g.set_color(chain_colors[IMP.atom.Chain(d1).get_id()])
107  g.set_name(get_nice_name(h) + "_skel")
108  w.add_geometry(g)
109 
110 
111 IMP.set_log_level(IMP.TERSE)
112 m = IMP.Model()
113 
114 # change to 46 to display all of them
115 hs = read(m, 3)
116 
117 # used to test of two molecules are touching one another
118 ps = IMP.core.KClosePairsPairScore(
119  IMP.core.SphereDistancePairScore(IMP.core.HarmonicUpperBound(10, 1)),
120  IMP.core.LeavesRefiner(IMP.atom.Hierarchy.get_traits()))
121 ps.set_log_level(IMP.SILENT)
122 
123 
124 print("creating rigid bodies")
125 base_chains = {}
126 for hc in IMP.atom.get_by_type(hs[0], IMP.atom.CHAIN_TYPE):
127  c = IMP.atom.Chain(hc)
128  base_chains[c.get_id()] = c
129 
130 for i, h in enumerate(hs):
131  for hc in IMP.atom.get_by_type(h, IMP.atom.CHAIN_TYPE):
132  c = IMP.atom.Chain(hc)
133  if h == hs[0]:
134  crb = IMP.atom.create_rigid_body(hc)
135  else:
136  # Make sure the rigid bodies have equivalent defining reference
137  # frames. If we just used IMP.atom.create_rigid_body, globular
138  # proteins are likely to have different axes computed when
139  # starting in different orientations
140  crb = IMP.atom.create_compatible_rigid_body(
141  hc, base_chains[c.get_id()])
142  print(" ", i)
143 
144 chains = IMP.atom.get_by_type(hs[0], IMP.atom.CHAIN_TYPE)
145 chains.sort(key=lambda x: IMP.core.XYZ(x).get_x() + IMP.core.XYZ(x).get_y())
146 chain_colors = {}
147 for i, c in enumerate(chains):
148  id = IMP.atom.Chain(c).get_id()
149  color = IMP.display.get_display_color(i)
150  # IMP.display.get_jet_color(f)
151  chain_colors[id] = color
152 
153 w = IMP.display.PymolWriter("markers.pym")
154 add_markers(hs[0], IMP.display.Color(1, 1, 1), w)
155 hso = hs[1:]
156 
157 
158 # sort them spatially so the colors are nicely arranged and allow one to
159 # visually connect the position of one end with that of the other
160 hso.sort(key=lambda h: IMP.core.XYZ(IMP.atom.Selection(
161  h, chain='I', residue_index=237).get_selected_particles()[0]).get_z())
162 print("adding markers", end=' ')
163 for i, h in enumerate(hso):
164  c = IMP.display.get_interpolated_rgb(
165  IMP.display.Color(1, 0, 0), IMP.display.Color(0, 0, 1), i / 50.)
166  add_markers(h, c, w)
167  print(" ", i)
168 w = IMP.display.PymolWriter("axis.pym")
169 print("adding axis", end=' ')
170 add_axis(hs[0], IMP.display.Color(1, 1, 1), w, chain_colors)
171 for i, h in enumerate(hs[1:]):
172  add_axis(h, None, w, chain_colors)
173  print(i, end=' ')
174 
175 w = IMP.display.PymolWriter("skeletons.pym")
176 add_skeleton(hs[0], IMP.display.Color(1, 1, 1), 5, w, chain_colors)
177 print("adding skeleton", end=' ')
178 for i, h in enumerate(hs[1:]):
179  add_skeleton(h, None, 1, w, chain_colors)
180  print(" ", i)