IMP  2.1.0
The Integrative Modeling Platform
kernel/nup84.py

Modify solve for a Nup84-like structure using a mix of rigid bodies and coarse grained models using crosslinking data. In addition, show how to visualize restraints and visualize the rejected conformations. Both are useful things to do when trying to figure out why optimization is not converging.

1 ## \example kernel/nup84.py
2 # Modify solve for a Nup84-like structure using a mix of rigid bodies
3 # and coarse grained models using crosslinking data. In
4 # addition, show how to visualize restraints and visualize the
5 # rejected conformations. Both are useful things to do when trying to
6 # figure out why optimization is not converging.
7 
8 import IMP
9 import IMP.atom
10 import IMP.container
11 import IMP.display
12 import IMP.statistics
13 import IMP.example
14 import os
15 import sys
16 
17 # not finished
18 IMP.base.add_bool_flag("run", "Whether to run the program")
19 
20 # parse command line arguments so, eg profiling can be used
21 IMP.base.setup_from_argv(sys.argv, "Nup84 example")
22 
23 if IMP.base.get_bool_flag("run") != "yes":
24  exit(0)
25 
26 # First we define some basic parameters for the modeling effort
27 
28 # the spring constant to use, it doesn't really matter
29 k = 10
30 # the target resolution for the representation
31 resolution = 100
32 # the box to perform everything in, make it flat as it is a 2D structure
34  IMP.algebra.Vector3D(300, 300, 50))
35 
36 # how many times to try to find a good solution
37 number_of_sampling_attempts = 1
38 number_of_mc_steps = 10000
39 
40 # Create a coarse grained protein with a given name, adding it to universe
41 
42 
43 def add_protein_from_length(model, name, residues, parent, restraints,
44  excluded_volume_particles, optimized_particles):
45  # Create a coarse grained protein with the passed residue information
46  h = IMP.atom.create_protein(model, name, resolution, residues)
47 
48  parent.add_child(h)
49  # Here, each of the constituent particles will be optimized independently
50  leaves = IMP.atom.get_leaves(h)
51  optimized_particles.extend(leaves)
52  excluded_volume_particles.extend(leaves)
53 
54  # Ensure that the various particles of the protein stay connected
56  for c in h.get_children()], k,
57  "connect " + name)
58 
59  if r:
60  # make sure there is an actual restraint
61  restraints.append(r)
62  r.set_maximum_score(k)
63 
64 # Create protein as a rigid body from a pdb file
65 
66 
67 def add_protein_from_pdb(model, name, file, parent, restraints,
68  excluded_volume_particles, optimized_particles):
69  # we should keep the original particles around so they get written
70 
71  # create an atomic representation from the pdb file
73  IMP.get_example_path(os.path.join("data", file)), model,
75  # extract the chain from the file
76  c = IMP.atom.Chain(IMP.atom.get_by_type(t, IMP.atom.CHAIN_TYPE)[0])
77  # there is no reason to use all atoms, just approximate the pdb shape
78  # instead
79  s = IMP.atom.create_simplified_along_backbone(c, resolution / 2.0, True)
80  s.set_name(name)
81  # tear down what is left
83  # make the simplified structure rigid
85  rb.set_coordinates_are_optimized(True)
86  optimized_particles.append(rb)
87  excluded_volume_particles.extend(s.get_children())
88  parent.add_child(s)
89 
90 
91 # Create protein as a rigid body from several pdb file
92 def add_protein_from_pdbs(model, name, files, parent, restraints,
93  excluded_volume_particles, optimized_particles):
95  for i, f in enumerate(files):
96  add_protein_from_pdb(model, name + str(i), f, h, restraints,
97  excluded_volume_particles, optimized_particles)
98  r = IMP.atom.create_connectivity_restraint([IMP.atom.Selection(c, hierarchy_types=[IMP.atom.FRAGMENT_TYPE])
99  for c in h.get_children()],
100  k, "connect " + name)
101  if r:
102  restraints.append(r)
103  r.set_maximum_score(k)
104 
105 # Create all the needed representation using the above functions
106 
107 
108 def create_representation(model):
109  restraints = []
110  optimized_particles = []
111  excluded_volume_particles = []
113  IMP.kernel.Particle(model, "the universe"))
114 
115  add_protein_from_length(model, "Nup85", 570, universe, restraints,
116  excluded_volume_particles, optimized_particles)
117 
118  # pin the c-terminus
119  ct = IMP.atom.Selection(universe, molecule="Nup85",
120  hierarchy_types=[IMP.atom.FRAGMENT_TYPE],
121  terminus=IMP.atom.Selection.C)
122  d = IMP.core.XYZ(ct.get_selected_particles()[0])
123  d.set_coordinates(IMP.algebra.Vector3D(0, 0, 0))
124  d.set_coordinates_are_optimized(False)
125 
126  add_protein_from_length(model, "Nup84", 460, universe, restraints,
127  excluded_volume_particles, optimized_particles)
128  add_protein_from_length(model, "Nup145C", 442, universe, restraints,
129  excluded_volume_particles, optimized_particles)
130  add_protein_from_length(
131  model, "Nup120", [0, 500, 761], universe, restraints,
132  excluded_volume_particles, optimized_particles)
133  add_protein_from_length(
134  model, "Nup133", [0, 450, 778, 1160], universe, restraints,
135  excluded_volume_particles, optimized_particles)
136  add_protein_from_pdb(model, "Seh1", "seh1.pdb", universe, restraints,
137  excluded_volume_particles, optimized_particles)
138  add_protein_from_pdb(model, "Sec13", "sec13.pdb", universe, restraints,
139  excluded_volume_particles, optimized_particles)
140  return universe, restraints, excluded_volume_particles, optimized_particles
141 
142 
143 def add_distance_restraint(selection0, selection1, name, restraints):
144  r = IMP.atom.create_distance_restraint(selection0, selection1, 0, k, name)
145  r.set_maximum_score(k)
146  restraints.append(r)
147 
148 
149 def encode_data_as_restraints(universe, restraints):
150  s0 = IMP.atom.Selection(
151  hierarchy=universe, hierarchy_types=[IMP.atom.FRAGMENT_TYPE],
152  molecule="Nup145C", residue_indexes=[(0, 423)])
153  s1 = IMP.atom.Selection(
154  hierarchy=universe, hierarchy_types=[IMP.atom.FRAGMENT_TYPE],
155  molecule="Nup84")
156  s2 = IMP.atom.Selection(
157  hierarchy=universe, hierarchy_types=[IMP.atom.FRAGMENT_TYPE],
158  molecule="Sec13")
160  [s0, s1, s2], k, "Nup145C Nup84 Sec13")
161  r.set_maximum_score(k)
162  restraints.append(r)
163 
164  add_distance_restraint(
165  IMP.atom.Selection(hierarchy=universe, molecule="Nup145C",
166  residue_indexes=[(0, 423)],
167  hierarchy_types=[
168  IMP.atom.FRAGMENT_TYPE]),
170  hierarchy=universe, molecule="Nup85",
171  hierarchy_types=[
172  IMP.atom.FRAGMENT_TYPE]),
173  "Num145C, Nup85", restraints)
174  add_distance_restraint(
175  IMP.atom.Selection(hierarchy=universe, molecule="Nup145C",
176  residue_indexes=[(0, 423)],
177  hierarchy_types=[
178  IMP.atom.FRAGMENT_TYPE]),
180  hierarchy=universe, molecule="Nup120",
181  residue_indexes=[(500, 762)],
182  hierarchy_types=[
183  IMP.atom.FRAGMENT_TYPE]),
184  "Nup145C Nup120", restraints)
185  add_distance_restraint(
186  IMP.atom.Selection(hierarchy=universe, molecule="Nup84",
187  hierarchy_types=[
188  IMP.atom.FRAGMENT_TYPE]),
190  hierarchy=universe, molecule="Nup133",
191  residue_indexes=[(778, 1160)],
192  hierarchy_types=[
193  IMP.atom.FRAGMENT_TYPE]),
194  "Nup84 Nup133", restraints)
195  add_distance_restraint(
196  IMP.atom.Selection(hierarchy=universe, molecule="Nup85",
197  hierarchy_types=[
198  IMP.atom.FRAGMENT_TYPE]),
200  hierarchy=universe, molecule="Seh1",
201  hierarchy_types=[
202  IMP.atom.FRAGMENT_TYPE]),
203  "Nup85 Seh1", restraints)
204  add_distance_restraint(
205  IMP.atom.Selection(hierarchy=universe, molecule="Nup145C",
206  residue_indexes=[(0, 423)],
207  hierarchy_types=[
208  IMP.atom.FRAGMENT_TYPE]),
210  hierarchy=universe, molecule="Sec13",
211  hierarchy_types=[
212  IMP.atom.FRAGMENT_TYPE]),
213  "Nup145C Sec13", restraints)
214 
215 
216 # find acceptable conformations of the model
217 def get_configurations(model, restraints, excluded_volume_particles, optimized_particles):
218  #cpc= IMP.container.ClosePairContainer(representation.get_particles(), 0, 10)
219  # evr= IMP.container.PairRestraint(IMP.core.SoftSpherePairScore(k), cpc,
220  # "Excluded Volume")
221  scale = .5
222  mc = IMP.core.MonteCarlo(model)
223  movers = []
224  for p in optimized_particles:
226  mover = IMP.core.RigidBodyMover(
227  p, IMP.core.XYZR(p).get_radius() * scale,
228  .2 * scale)
229  movers.append(mover)
230  else:
231  mover = IMP.core.BallMover(
232  [p], IMP.core.XYZR(p).get_radius() * scale)
233  movers.append(mover)
234  serial_mover = IMP.core.SerialMover(movers)
235  mc.add_mover(serial_mover)
236  scoring_function = IMP.core.IncrementalScoringFunction(
237  optimized_particles, restraints)
238  scoring_function.add_close_pair_score(IMP.core.SoftSpherePairScore(k), 0.0,
239  excluded_volume_particles)
240 
241  configuration_set = IMP.ConfigurationSet(model)
242  # must write our own sampler as IMP.core.MCCGSampler doesn't handle rigid
243  # bodies
244  for i in range(number_of_sampling_attempts):
245  for p in optimized_particles:
246  IMP.core.XYZ(p).set_coordinates(
248  mc.optimize(number_of_mc_steps)
249  if scoring_function.get_had_good_score():
250  configuration_set.save()
251  return configuration_set
252 
253 
254 model = IMP.kernel.Model()
255 universe, restraints, excluded_volume_particles, optimized_particles = create_representation(
256  model)
257 encode_data_as_restraints(universe, restraints)
258 
259 configuration_set = get_configurations(model, restraints,
260  excluded_volume_particles,
261  optimized_particles)
262 
263 print "Found", configuration_set.get_number_of_configurations(), "good configurations"
264 
265 # now lets save them all to a file
266 rmf_file_name = IMP.base.create_temporary_file_name("nup84", ".rmf")
267 rmf = RMF.create_rmf_file(rmf_file_name)
268 
269 # we want to see the scores of the various restraints also
270 IMP.rmf.add_restraints(rmf, restraints)
271 # and the actual structures
272 IMP.rmf.add_hierarchy(rmf, universe)
273 
274 for i in range(0, configuration_set.get_number_of_configurations()):
275  configuration_set.load_configuration(i)
276  # align all the configurations with the first so they display nicely
277  # if we want to be fancy we can account for flips too
278  if i == 0:
279  base_coords = [IMP.core.XYZ(p).get_coordinates()
280  for p in optimized_particles]
281  else:
282  tr = IMP.algebra.get_transform_taking_first_to_second(
283  optimized_particles, base_coords)
284  IMP.core.transform(optimized_particles, tr)
285  # update the restraint scores
286  sf.evaluate(False)
287  IMP.rmf.save_frame(rmf, i)
288 
289 print "You can now open", rmf_file_name, "in chimera"