IMP logo
IMP Reference Guide  2.8.0
The Integrative Modeling Platform
samplers.py
1 """@namespace IMP.pmi.samplers
2  Sampling of the system.
3 """
4 
5 from __future__ import print_function
6 import IMP
7 import IMP.core
8 from IMP.pmi.tools import get_restraint_set
9 
10 class _SerialReplicaExchange(object):
11  """Dummy replica exchange class used in non-MPI builds.
12  It should act similarly to IMP.mpi.ReplicaExchange on a single processor.
13  """
14  def __init__(self):
15  self.__params = {}
16  def get_number_of_replicas(self):
17  return 1
18  def create_temperatures(self, tmin, tmax, nrep):
19  return [tmin]
20  def get_my_index(self):
21  return 0
22  def set_my_parameter(self, key, val):
23  self.__params[key] = val
24  def get_my_parameter(self, key):
25  return self.__params[key]
26  def get_friend_index(self, step):
27  return 0
28  def get_friend_parameter(self, key, findex):
29  return self.get_my_parameter(key)
30  def do_exchange(self, myscore, fscore, findex):
31  return False
32  def set_was_used(self,was_used):
33  self.was_used = was_used
34 
35 
36 class MonteCarlo(object):
37  """Sample using Monte Carlo"""
38 
39  # check that isd is installed
40  try:
41  import IMP.isd
42  isd_available = True
43  except ImportError:
44  isd_available = False
45 
46  def __init__(self, m, objects=None, temp=1.0, filterbyname=None):
47  """Setup Monte Carlo sampling
48  @param m The IMP Model
49  @param objects What to sample. Use flat list of particles or
50  (deprecated) 'MC Sample Objects' from PMI1
51  @param temp The MC temperature
52  @param filterbyname Not used
53  """
54  self.losp = [
55  "Rigid_Bodies",
56  "Floppy_Bodies",
57  "Nuisances",
58  "X_coord",
59  "Weights"
60  "Surfaces"]
61  self.simulated_annealing = False
62  self.selfadaptive = False
63  # that is -1 because mc has not yet run
64  self.nframe = -1
65  self.temp = temp
66  self.mvs = []
67  self.mvslabels = []
68  self.label = "None"
69  self.m = m
70  self.movers_data={}
71 
72  # check if using PMI1 or just passed a list of movers
73  gather_objects = False
74  try:
75  objects[0].get_particles_to_sample()
76  gather_objects = True
77  except:
78  self.mvs = objects
79 
80  if gather_objects:
81  for ob in objects:
82  pts = ob.get_particles_to_sample()
83  for k in pts.keys():
84 
85  if "Rigid_Bodies" in k:
86  mvs = self.get_rigid_body_movers(
87  pts[k][0],
88  pts[k][1],
89  pts[k][2])
90  for mv in mvs:
91  mv.set_name(k)
92  self.mvs += mvs
93 
94  if "SR_Bodies" in k:
95  print(len(pts[k]))
96  mvs = self.get_super_rigid_body_movers(
97  pts[k][0],
98  pts[k][1],
99  pts[k][2])
100  for mv in mvs:
101  mv.set_name(k)
102  self.mvs += mvs
103 
104  if "Floppy_Bodies" in k:
105  mvs = self.get_floppy_body_movers(pts[k][0], pts[k][1])
106  for mv in mvs:
107  mv.set_name(k)
108  self.mvs += mvs
109 
110  if "X_coord" in k:
111  mvs = self.get_X_movers(pts[k][0], pts[k][1])
112  for mv in mvs:
113  mv.set_name(k)
114  self.mvs += mvs
115 
116  if "Nuisances" in k:
117  if not self.isd_available:
118  raise ValueError("isd module needed to use nuisances")
119  mvs = self.get_nuisance_movers(pts[k][0], pts[k][1])
120  for mv in mvs:
121  mv.set_name(k)
122  self.mvs += mvs
123 
124  if "Weights" in k:
125  if not self.isd_available:
126  raise ValueError("isd module needed to use weights")
127  mvs = self.get_weight_movers(pts[k][0], pts[k][1])
128  for mv in mvs:
129  mv.set_name(k)
130  self.mvs += mvs
131 
132  if "Surfaces" in k:
133  mvs = self.get_surface_movers(
134  pts[k][0],
135  pts[k][1],
136  pts[k][2],
137  pts[k][3])
138  for mv in mvs:
139  mv.set_name(k)
140  self.mvs += mvs
141 
142  # SerialMover
143  self.smv = IMP.core.SerialMover(self.mvs)
144 
145  self.mc = IMP.core.MonteCarlo(self.m)
146  self.mc.set_scoring_function(get_restraint_set(self.m))
147  self.mc.set_return_best(False)
148  self.mc.set_kt(self.temp)
149  self.mc.add_mover(self.smv)
150 
151  def set_kt(self, temp):
152  self.temp = temp
153  self.mc.set_kt(temp)
154 
155  def get_mc(self):
156  return self.mc
157 
158  def set_scoring_function(self, objectlist):
159  rs = IMP.RestraintSet(self.m, 1.0, 'sfo')
160  for ob in objectlist:
161  rs.add_restraint(ob.get_restraint())
163  self.mc.set_scoring_function(sf)
164 
165  def set_simulated_annealing(
166  self,
167  min_temp,
168  max_temp,
169  min_temp_time,
170  max_temp_time):
171  self.simulated_annealing = True
172  self.tempmin = min_temp
173  self.tempmax = max_temp
174  self.timemin = min_temp_time
175  self.timemax = max_temp_time
176 
177  def set_self_adaptive(self, isselfadaptive=True):
178  self.selfadaptive = isselfadaptive
179 
181  '''
182  Return a dictionary with the mover parameters for nuisance parameters
183  '''
184  output = {}
185  for i in range(self.get_number_of_movers()):
186  mv = self.smv.get_mover(i)
187  name = mv.get_name()
188  if "Nuisances" in name:
189  stepsize = IMP.core.NormalMover.get_from(mv).get_sigma()
190  output[name] = stepsize
191  return output
192 
193  def get_number_of_movers(self):
194  return len(self.smv.get_movers())
195 
196  def get_particle_types():
197  return self.losp
198 
199  def optimize(self, nstep):
200  self.nframe += 1
201  self.mc.optimize(nstep * self.get_number_of_movers())
202 
203  # apply simulated annealing protocol
204  if self.simulated_annealing:
205  self.temp = self.temp_simulated_annealing()
206  self.mc.set_kt(self.temp)
207 
208  # apply self adaptive protocol
209  if self.selfadaptive:
210  for i, mv in enumerate(self.mvs):
211 
212  mvacc = mv.get_number_of_accepted()
213  mvprp = mv.get_number_of_proposed()
214  if mv not in self.movers_data:
215  accept = float(mvacc) / float(mvprp)
216  self.movers_data[mv]=(mvacc,mvprp)
217  else:
218  oldmvacc,oldmvprp=self.movers_data[mv]
219  accept = float(mvacc-oldmvacc) / float(mvprp-oldmvprp)
220  self.movers_data[mv]=(mvacc,mvprp)
221  if accept < 0.05: accept = 0.05
222  if accept > 1.0: accept = 1.0
223 
224  if type(mv) is IMP.core.NormalMover:
225  stepsize = mv.get_sigma()
226  if 0.4 > accept or accept > 0.6:
227  mv.set_sigma(stepsize * 2 * accept)
228 
229  if type(mv) is IMP.isd.WeightMover:
230  stepsize = mv.get_radius()
231  if 0.4 > accept or accept > 0.6:
232  mv.set_radius(stepsize * 2 * accept)
233 
234  if type(mv) is IMP.core.RigidBodyMover:
235  mr=mv.get_maximum_rotation()
236  mt=mv.get_maximum_translation()
237  if 0.4 > accept or accept > 0.6:
238  mv.set_maximum_rotation(mr * 2 * accept)
239  mv.set_maximum_translation(mt * 2 * accept)
240 
241  if type(mv) is IMP.pmi.TransformMover:
242  mr=mv.get_maximum_rotation()
243  mt=mv.get_maximum_translation()
244  if 0.4 > accept or accept > 0.6:
245  mv.set_maximum_rotation(mr * 2 * accept)
246  mv.set_maximum_translation(mt * 2 * accept)
247 
248  if type(mv) is IMP.core.BallMover:
249  mr=mv.get_radius()
250  if 0.4 > accept or accept > 0.6:
251  mv.set_radius(mr * 2 * accept)
252 
253  @IMP.deprecated_method("2.5", "Use optimize() instead.")
254  def run(self, *args, **kwargs):
255  self.optimize(*args, **kwargs)
256 
257  def get_nuisance_movers(self, nuisances, maxstep):
258  mvs = []
259  for nuisance in nuisances:
260  print(nuisance, maxstep)
261  mvs.append(
262  IMP.core.NormalMover([nuisance],
263  IMP.FloatKeys([IMP.FloatKey("nuisance")]),
264  maxstep))
265  return mvs
266 
267  def get_rigid_body_movers(self, rbs, maxtrans, maxrot):
268  mvs = []
269  for rb in rbs:
270  mvs.append(IMP.core.RigidBodyMover(rb.get_model(), rb,
271  maxtrans, maxrot))
272  return mvs
273 
274  def get_super_rigid_body_movers(self, rbs, maxtrans, maxrot):
275  mvs = []
276  for rb in rbs:
277  if len(rb) == 2:
278  # normal Super Rigid Body
279  srbm = IMP.pmi.TransformMover(self.m, maxtrans, maxrot)
280  elif len(rb) == 3:
281  if type(rb[2]) == tuple and type(rb[2][0]) == float \
282  and type(rb[2][1]) == float and type(rb[2][2]) == float \
283  and len(rb[2])== 3:
284  # super rigid body with 2D rotation, rb[2] is the axis
285  srbm = IMP.pmi.TransformMover(
286  self.m,
287  IMP.algebra.Vector3D(rb[2]),
288  maxtrans,
289  maxrot)
290  #elif type(rb[2]) == tuple and type(rb[2][0]) == IMP.Particle \
291  # and type(rb[2][1]) == IMP.Particle and len(rb[2])== 2:
292  # # super rigid body with bond rotation
293 
294  # srbm = IMP.pmi.TransformMover(
295  # self.m,
296  # rb[2][0],rb[2][1],
297  # 0, #no translation
298  # maxrot)
299  else:
300  print("Setting up a super rigid body with wrong parameters")
301  raise
302 
303  for xyz in rb[0]:
304  srbm.add_xyz_particle(xyz)
305  for rb in rb[1]:
306  srbm.add_rigid_body_particle(rb)
307  mvs.append(srbm)
308  return mvs
309 
310  def get_floppy_body_movers(self, fbs, maxtrans):
311  mvs = []
312  for fb in fbs:
313  # check is that is a rigid body member:
315  # if so force the particles to move anyway
316  floatkeys = [IMP.FloatKey(4), IMP.FloatKey(5), IMP.FloatKey(6)]
317  for fk in floatkeys:
318  fb.set_is_optimized(fk, True)
319  mvs.append(
320  IMP.core.BallMover(fb.get_model(), fb,
321  IMP.FloatKeys(floatkeys),
322  maxtrans))
323  else:
324  # otherwise use the normal ball mover
325  mvs.append(IMP.core.BallMover(fb.get_model(), fb, maxtrans))
326  return mvs
327 
328  def get_X_movers(self, fbs, maxtrans):
329  mvs = []
330  Xfloatkey = IMP.core.XYZ.get_xyz_keys()[0]
331  for fb in fbs:
332  # check is that is a rigid body member:
334  raise ValueError("particle is part of a rigid body")
335  else:
336  # otherwise use the normal ball mover
337  mvs.append(IMP.core.NormalMover([fb], [Xfloatkey], maxtrans))
338  return mvs
339 
340  def get_weight_movers(self, weights, maxstep):
341  mvs = []
342  for weight in weights:
343  if(weight.get_number_of_states() > 1):
344  mvs.append(IMP.isd.WeightMover(weight, maxstep))
345  return mvs
346 
347  def get_surface_movers(self, surfaces, maxtrans, maxrot, refprob):
348  mvs = []
349  for surface in surfaces:
350  mvs.append(IMP.core.SurfaceMover(surface, maxtrans, maxrot,
351  refprob))
352  return mvs
353 
354  def temp_simulated_annealing(self):
355  if self.nframe % (self.timemin + self.timemax) < self.timemin:
356  value = 0.0
357  else:
358  value = 1.0
359  temp = self.tempmin + (self.tempmax - self.tempmin) * value
360  return temp
361 
362  def set_label(self, label):
363  self.label = label
364 
365  def get_frame_number(self):
366  return self.nframe
367 
368  def get_output(self):
369  output = {}
370  acceptances = []
371  for i, mv in enumerate(self.smv.get_movers()):
372  mvname = mv.get_name()
373  mvacc = mv.get_number_of_accepted()
374  mvprp = mv.get_number_of_proposed()
375  try:
376  mvacr = float(mvacc) / float(mvprp)
377  except:
378  mvacr = 0.0
379  output["MonteCarlo_Acceptance_" +
380  mvname + "_" + str(i)] = str(mvacr)
381  if "Nuisances" in mvname:
382  output["MonteCarlo_StepSize_" + mvname + "_" +
383  str(i)] = str(IMP.core.NormalMover.get_from(mv).get_sigma())
384  if "Weights" in mvname:
385  output["MonteCarlo_StepSize_" + mvname + "_" +
386  str(i)] = str(IMP.isd.WeightMover.get_from(mv).get_radius())
387  output["MonteCarlo_Temperature"] = str(self.mc.get_kt())
388  output["MonteCarlo_Nframe"] = str(self.nframe)
389  return output
390 
391 
392 class MolecularDynamics(object):
393  """Sample using molecular dynamics"""
394 
395  def __init__(self,m,objects,kt,gamma=0.01,maximum_time_step=1.0,sf=None):
396  """Setup MD
397  @param m The IMP Model
398  @param objects What to sample. Use flat list of particles or (deprecated) 'MD Sample Objects' from PMI1
399  @param kt Temperature
400  @param gamma Viscosity parameter
401  @param maximum_time_step MD max time step
402  """
403  self.m=m
404 
405  # check if using PMI1 objects dictionary, or just list of particles
406  try:
407  for obj in objects:
408  to_sample=obj.get_particles_to_sample()['Floppy_Bodies_SimplifiedModel'][0]
409  except:
410  to_sample = objects
411 
412  self.ltstate=IMP.atom.LangevinThermostatOptimizerState(self.m,to_sample,
413  kt/0.0019872041,
414  gamma)
415  self.md = IMP.atom.MolecularDynamics(self.m)
416  self.md.set_maximum_time_step(maximum_time_step)
417  if sf:
418  self.md.set_scoring_function(sf)
419  else:
420  self.md.set_scoring_function(get_restraint_set(self.m))
421  self.md.add_optimizer_state(self.ltstate)
422  self.simulated_annealing = False
423  self.nframe = -1
424  def set_kt(self,kt):
425  temp=kt/0.0019872041
426  self.ltstate.set_temperature(temp)
427  self.md.assign_velocities(temp)
428 
429  def set_simulated_annealing(self, min_temp, max_temp, min_temp_time,
430  max_temp_time):
431  self.simulated_annealing = True
432  self.tempmin = min_temp
433  self.tempmax = max_temp
434  self.timemin = min_temp_time
435  self.timemax = max_temp_time
436 
437  def temp_simulated_annealing(self):
438  if self.nframe % (self.timemin + self.timemax) < self.timemin:
439  value = 0.0
440  else:
441  value = 1.0
442  temp = self.tempmin + (self.tempmax - self.tempmin) * value
443  return temp
444 
445  def set_gamma(self,gamma):
446  self.ltstate.set_gamma(gamma)
447 
448  def optimize(self,nsteps):
449  # apply simulated annealing protocol
450  self.nframe+=1
451  if self.simulated_annealing:
452  self.temp = self.temp_simulated_annealing()
453  self.set_kt(self.temp)
454  self.md.optimize(nsteps)
455 
456  def get_output(self):
457  output={}
458  output["MolecularDynamics_KineticEnergy"]=str(self.md.get_kinetic_energy())
459  return output
460 
461 class ConjugateGradients(object):
462  """Sample using conjugate gradients"""
463 
464  def __init__(self, m, objects):
465  self.m = m
466  self.nframe = -1
467  self.cg = IMP.core.ConjugateGradients(self.m)
468  self.cg.set_scoring_function(get_restraint_set(self.m))
469 
470  def set_label(self, label):
471  self.label = label
472 
473  def get_frame_number(self):
474  return self.nframe
475 
476  @IMP.deprecated_method("2.5", "Use optimize() instead.")
477  def run(self, *args, **kwargs):
478  self.optimize(*args, **kwargs)
479 
480  def optimize(self, nstep):
481  self.nframe += 1
482  self.cg.optimize(nstep)
483 
484  def set_scoring_function(self, objectlist):
485  rs = IMP.RestraintSet(self.m, 1.0, 'sfo')
486  for ob in objectlist:
487  rs.add_restraint(ob.get_restraint())
489  self.cg.set_scoring_function(sf)
490 
491  def get_output(self):
492  output = {}
493  acceptances = []
494  output["ConjugatedGradients_Nframe"] = str(self.nframe)
495  return output
496 
497 
498 
499 
500 
501 
502 class ReplicaExchange(object):
503  """Sample using replica exchange"""
504 
505  def __init__(
506  self,
507  model,
508  tempmin,
509  tempmax,
510  samplerobjects,
511  test=True,
512  replica_exchange_object=None):
513  '''
514  samplerobjects can be a list of MonteCarlo or MolecularDynamics
515  '''
516 
517 
518  self.m = model
519  self.samplerobjects = samplerobjects
520  # min and max temperature
521  self.TEMPMIN_ = tempmin
522  self.TEMPMAX_ = tempmax
523 
524  if replica_exchange_object is None:
525  # initialize Replica Exchange class
526  try:
527  import IMP.mpi
528  print('ReplicaExchange: MPI was found. Using Parallel Replica Exchange')
529  self.rem = IMP.mpi.ReplicaExchange()
530  except ImportError:
531  print('ReplicaExchange: Could not find MPI. Using Serial Replica Exchange')
532  self.rem = _SerialReplicaExchange()
533 
534  else:
535  # get the replica exchange class instance from elsewhere
536  print('got existing rex object')
537  self.rem = replica_exchange_object
538 
539  # get number of replicas
540  nproc = self.rem.get_number_of_replicas()
541 
542  if nproc % 2 != 0 and test == False:
543  raise Exception("number of replicas has to be even. set test=True to run with odd number of replicas.")
544  # create array of temperatures, in geometric progression
545  temp = self.rem.create_temperatures(
546  self.TEMPMIN_,
547  self.TEMPMAX_,
548  nproc)
549  # get replica index
550  self.temperatures = temp
551 
552  myindex = self.rem.get_my_index()
553  # set initial value of the parameter (temperature) to exchange
554  self.rem.set_my_parameter("temp", [self.temperatures[myindex]])
555  for so in self.samplerobjects:
556  so.set_kt(self.temperatures[myindex])
557  self.nattempts = 0
558  self.nmintemp = 0
559  self.nmaxtemp = 0
560  self.nsuccess = 0
561 
562  def get_temperatures(self):
563  return self.temperatures
564 
565  def get_my_temp(self):
566  return self.rem.get_my_parameter("temp")[0]
567 
568  def get_my_index(self):
569  return self.rem.get_my_index()
570 
571  def swap_temp(self, nframe, score=None):
572  if score is None:
573  score = self.m.evaluate(False)
574  # get my replica index and temperature
575  myindex = self.rem.get_my_index()
576  mytemp = self.rem.get_my_parameter("temp")[0]
577 
578  if mytemp == self.TEMPMIN_:
579  self.nmintemp += 1
580 
581  if mytemp == self.TEMPMAX_:
582  self.nmaxtemp += 1
583 
584  # score divided by kbt
585  myscore = score / mytemp
586 
587  # get my friend index and temperature
588  findex = self.rem.get_friend_index(nframe)
589  ftemp = self.rem.get_friend_parameter("temp", findex)[0]
590  # score divided by kbt
591  fscore = score / ftemp
592 
593  # try exchange
594  flag = self.rem.do_exchange(myscore, fscore, findex)
595 
596  self.nattempts += 1
597  # if accepted, change temperature
598  if (flag):
599  for so in self.samplerobjects:
600  so.set_kt(ftemp)
601  self.nsuccess += 1
602 
603  def get_output(self):
604  output = {}
605  acceptances = []
606  if self.nattempts != 0:
607  output["ReplicaExchange_SwapSuccessRatio"] = str(
608  float(self.nsuccess) / self.nattempts)
609  output["ReplicaExchange_MinTempFrequency"] = str(
610  float(self.nmintemp) / self.nattempts)
611  output["ReplicaExchange_MaxTempFrequency"] = str(
612  float(self.nmaxtemp) / self.nattempts)
613  else:
614  output["ReplicaExchange_SwapSuccessRatio"] = str(0)
615  output["ReplicaExchange_MinTempFrequency"] = str(0)
616  output["ReplicaExchange_MaxTempFrequency"] = str(0)
617  output["ReplicaExchange_CurrentTemp"] = str(self.get_my_temp())
618  return output
619 
620 
621 
622 class MPI_values(object):
623  def __init__(self,replica_exchange_object=None):
624  """Query values (ie score, and others)
625  from a set of parallel jobs"""
626 
627  if replica_exchange_object is None:
628  # initialize Replica Exchange class
629  try:
630  import IMP.mpi
631  print('MPI_values: MPI was found. Using Parallel Replica Exchange')
632  self.rem = IMP.mpi.ReplicaExchange()
633  except ImportError:
634  print('MPI_values: Could not find MPI. Using Serial Replica Exchange')
635  self.rem = _SerialReplicaExchange()
636 
637  else:
638  # get the replica exchange class instance from elsewhere
639  print('got existing rex object')
640  self.rem = replica_exchange_object
641 
642  def set_value(self,name,value):
643  self.rem.set_my_parameter(name,[value])
644 
645  def get_values(self,name):
646  values=[]
647  for i in range(self.rem.get_number_of_replicas()):
648  v=self.rem.get_friend_parameter(name, i)[0]
649  values.append(v)
650  return values
651 
652  def get_percentile(self,name):
653  value=self.rem.get_my_parameter(name)[0]
654  values=sorted(self.get_values(name))
655  ind=values.index(value)
656  percentile=float(ind)/len(values)
657  return percentile
658 
659 
660 
661 class PyMCMover(object):
662  # only works if the sampled particles are rigid bodies
663 
664  def __init__(self, representation, mcchild, n_mc_steps):
665 
666  # mcchild must be pmi.samplers.MonteCarlo
667  # representation must be pmi.representation
668 
669  self.rbs = representation.get_rigid_bodies()
670 
671  self.mc = mcchild
672  self.n_mc_steps = n_mc_steps
673 
674  def store_move(self):
675  # get all xyz coordinates of all rigid bodies of all copies
676  self.oldcoords = []
677  for copy in self.rbs:
678  crd = []
679  for rb in copy:
680  crd.append(rb.get_reference_frame())
681  self.oldcoords.append(crd)
682 
683  def propose_move(self, prob):
684  self.mc.run(self.n_mc_steps)
685 
686  def reset_move(self):
687  # reset each copy to crd
688  for copy, crd in zip(self.rbs, self.oldcoords):
689  for rb, ref in zip(copy, crd):
690  rb.set_reference_frame(ref)
691 
692  def get_number_of_steps(self):
693  return self.n_mc_steps
694 
695  def set_number_of_steps(self, nsteps):
696  self.n_mc_steps = nsteps
697 
698 
699 class PyMC(object):
700 
701  def __init__(self, model):
702  from math import exp
703  import random
704 
705  self.m = model
706  self.restraints = None
707  self.first_call = True
708  self.nframe = -1
709 
710  def add_mover(self, mv):
711  self.mv = mv
712 
713  def set_kt(self, kT):
714  self.kT = kT
715 
716  def set_return_best(self, thing):
717  pass
718 
719  def set_move_probability(self, thing):
720  pass
721 
722  def get_energy(self):
723  if self.restraints:
724  pot = sum([r.evaluate(False) for r in self.restraints])
725  else:
726  pot = self.m.evaluate(False)
727  return pot
728 
729  def metropolis(self, old, new):
730  deltaE = new - old
731  print(": old %f new %f deltaE %f new_epot: %f" % (old, new, deltaE,
732  self.m.evaluate(
733  False)), end=' ')
734  kT = self.kT
735  if deltaE < 0:
736  return True
737  else:
738  return exp(-deltaE / kT) > random.uniform(0, 1)
739 
740  def optimize(self, nsteps):
741  self.naccept = 0
742  self.nframe += 1
743  print("=== new MC call")
744  # store initial coordinates
745  if self.first_call:
746  self.mv.store_move()
747  self.first_call = False
748  for i in range(nsteps):
749  print("MC step %d " % i, end=' ')
750  # get total energy
751  old = self.get_energy()
752  # make a MD move
753  self.mv.propose_move(1)
754  # get new total energy
755  new = self.get_energy()
756  if self.metropolis(old, new):
757  # move was accepted: store new conformation
758  self.mv.store_move()
759  self.naccept += 1
760  print("accepted ")
761  else:
762  # move rejected: restore old conformation
763  self.mv.reset_move()
764  print(" ")
765 
766  def get_number_of_forward_steps(self):
767  return self.naccept
768 
769  def set_restraints(self, restraints):
770  self.restraints = restraints
771 
772  def set_scoring_function(self, objects):
773  # objects should be pmi.restraints
774  rs = IMP.RestraintSet(self.m, 1.0, 'sfo')
775  for ob in objects:
776  rs.add_restraint(ob.get_restraint())
777  self.set_restraints([rs])
778 
779  def get_output(self):
780  output = {}
781  acceptances = []
782  output["PyMC_Temperature"] = str(self.kT)
783  output["PyMC_Nframe"] = str(self.nframe)
784  return output
785 
786 
787 # 3
def __init__
samplerobjects can be a list of MonteCarlo or MolecularDynamics
Definition: samplers.py:505
A Monte Carlo optimizer.
Definition: MonteCarlo.h:45
A class to implement Hamiltonian Replica Exchange.
Maintains temperature during molecular dynamics.
Sample using molecular dynamics.
Definition: samplers.py:392
def get_nuisance_movers_parameters
Return a dictionary with the mover parameters for nuisance parameters.
Definition: samplers.py:180
Miscellaneous utilities.
Definition: tools.py:1
Modify the transformation of a rigid body.
Simple conjugate gradients optimizer.
Sample using conjugate gradients.
Definition: samplers.py:461
Move continuous particle variables by perturbing them within a ball.
Modify a surface orientation.
Definition: SurfaceMover.h:34
Object used to hold a set of restraints.
Definition: RestraintSet.h:36
Simple molecular dynamics optimizer.
def deprecated_method
Python decorator to mark a method as deprecated.
Definition: __init__.py:9599
Code that uses the MPI parallel library.
A mover that perturbs a Weight particle.
Definition: WeightMover.h:20
Modify the transformation of a rigid body.
def __init__
Setup Monte Carlo sampling.
Definition: samplers.py:46
Modify a set of continuous variables using a normal distribution.
Definition: NormalMover.h:20
Basic functionality that is expected to be used by a wide variety of IMP users.
Sample using Monte Carlo.
Definition: samplers.py:36
The general base class for IMP exceptions.
Definition: exception.h:49
VectorD< 3 > Vector3D
Definition: VectorD.h:395
static const FloatKeys & get_xyz_keys()
Get a vector containing the keys for x,y,z.
Applies a list of movers one at a time.
Definition: SerialMover.h:23
static bool get_is_setup(const IMP::ParticleAdaptor &p)
Definition: rigid_bodies.h:640
Sample using replica exchange.
Definition: samplers.py:502
Inferential scoring building on methods developed as part of the Inferential Structure Determination ...