2 from __future__
import print_function
18 from math
import sqrt, exp
19 from io
import StringIO
21 from random
import random
23 kB = (1.381 * 6.02214) / 4184.0
28 def __init__(self, cont, md, n_md_steps):
31 self.n_md_steps = n_md_steps
38 for i
in range(self.cont.get_number_of_particles()):
39 p = self.cont.get_particle(i)
42 self.oldcoords.append(d.get_coordinates())
44 def propose_move(self, prob):
45 self.md.optimize(self.n_md_steps)
48 for i
in range(self.cont.get_number_of_particles()):
49 p = self.cont.get_particle(i)
53 d.set_coordinates(self.oldcoords[i])
55 def get_number_of_steps(self):
56 return self.n_md_steps
58 def set_number_of_steps(self, nsteps):
59 self.n_md_steps = nsteps
65 def __init__(self, model):
68 def add_mover(self, mv):
74 def set_return_best(self, thing):
77 def set_move_probability(self, thing):
81 pot = self.m.evaluate(
False)
83 kin = self.mv.md.get_kinetic_energy()
86 def metropolis(self, old, new):
89 print(
": old %f new %f deltaE %f new_epot: %f" % (old, new, deltaE,
96 return exp(-deltaE / kT) > random()
98 def optimize(self, nsteps):
100 print(
"=== new MC call")
103 for i
in range(nsteps):
104 print(
"MC step %d " % i, end=
' ')
106 self.mv.md.assign_velocities(self.kT / kB)
108 old = self.get_energy()
110 self.mv.propose_move(1)
112 new = self.get_energy()
113 if self.metropolis(old, new):
123 def get_number_of_forward_steps(self):
129 """nonspecific methods used across all shared function objects.
131 - Their name starts with the name of the parent function (e.g.
133 - they don't store anything in the class, but instead
134 return all created objects.
135 Exceptions: the model, which is self._m
136 the statistics class, which is self.stat
137 - they store what they have to store in the model (e.g. restraints)
138 - they don't print out anything except for long routines (e.g. NOE
140 - the prior RestraintSet is added to the model when it is created, so
141 that when it is passed to another function, it is not added twice.
147 def set_checklevel(self, value):
150 def set_loglevel(self, value):
153 def m(self, name, *args, **kw):
154 "wrapper to call methods of m"
155 func = getattr(self._m, name)
156 return func(*args, **kw)
159 "moves to wd and creates model"
165 self, initpdb, top, par, selector, pairscore,
166 ff_temp=300.0, disulfides=
None, representation=
'custom'):
167 """creates a CHARMM protein representation.
168 creates the charmm force field, bonded and nonbonded.
169 - initpdb: initial structure in pdb format
170 - top is a CHARMM top.lib, read if representation=='custom' (default)
171 - par is a CHARMM par.lib
172 - selector is an instance of
173 one of the selectors of IMP.atom, for example
174 IMP.atom.NonWaterNonHydrogenPDBSelector().
175 - pairscore is an instance of a Pair Score to score the interaction
176 between two atoms. usually, it's either
177 LennardJonesDistancePairScore(0,1) or
178 RepulsiveDistancePairScore(0,1)
179 - ff_temp is the temperature at which the force field should be
181 - disulfides: if not None, a list of tuples corresponding to residue
182 numbers that should be cross-linked. Residues should be
183 cysteines, and residue numbering should start at 0.
184 - representation: 'full' : all-atom CHARMM force field
185 'heavy': heavy atom CHARMM forcefield with polar H
186 'calpha': only C alphas, ball and stick model with
187 bondlength 3.78 angstrom, beads at VdW
188 contact, and harmonic restraint between
190 'custom' (default) : read given CHARMM top and par
192 Returns: prot, ff, rsb, rs
194 - ff: the force field
195 - rsb: the RestraintSet on bonded interactions
196 - rs: the RestraintSet on nonbonded interactions. Both are weighted
201 if not prot.get_is_valid(
True):
202 raise ValueError(
"invalid hierarchy!")
203 if representation ==
'custom':
206 elif representation ==
'heavy':
208 elif representation ==
'full':
210 elif representation ==
'calpha':
213 raise NotImplementedError(representation)
214 if representation ==
'calpha':
215 print(
"setting up simplified C alpha force field")
217 for ca
in IMP.atom.get_by_type(prot, IMP.atom.ATOM_TYPE):
222 for chain
in prot.get_children():
223 residues = [(chain.get_child(i), chain.get_child(i + 1))
224 for i
in range(chain.get_number_of_children() - 1)]
225 residues = [(i.get_child(0).get_particle(),
226 j.get_child(0).get_particle())
227 for (i, j)
in residues]
228 pairs.extend(residues)
248 rsb.add_restraint(br)
249 rsb.set_weight(1.0 / (kB * ff_temp))
251 nonbonded_pair_filter.set_bonds(bonds)
254 print(
"setting up CHARMM forcefield")
258 topology = ff.create_topology(prot)
262 topology.apply_default_patches()
265 s = topology.get_segment(0)
266 dis = ff.get_patch(
'DISU')
267 for (i, j)
in disulfides:
270 r0 = s.get_residue(i)
271 r1 = s.get_residue(j)
273 r0.set_patched(
False)
275 r1.set_patched(
False)
277 print(
"added disulfide bridge between cysteines %d and %d" % (i, j))
283 topology.setup_hierarchy(prot)
295 ff.add_well_depths(prot)
296 nonbonded_pair_filter = r.get_pair_filter()
298 atoms = IMP.atom.get_by_type(prot, IMP.atom.ATOM_TYPE)
309 nbl.add_pair_filter(nonbonded_pair_filter)
314 return prot, ff, rsb, rs
317 """sets up a Scale particle to the initial default value. It can
318 optionally be constrained between two positive bounds, or else its
319 range is 0 to infinity.
324 scale.set_lower(lower)
326 scale.set_upper(upper)
330 """given a list of scales, returns a RestraintSet('prior') with weight
331 1.0 that contains a list of vonMisesKappaJeffreysRestraint on each scale.
332 If argument prior_rs is used, add them to that RestraintSet instead.
336 prior_rs.set_weight(1.0)
342 """given a list of scales, returns a RestraintSet('prior') with weight
343 1.0 that contains a list of JeffreysRestraint on each scale.
344 If argument prior_rs is used, add them to that RestraintSet instead.
350 prior_rs.add_restraint(jr)
354 """given a list of scales, returns a RestraintSet('prior') with weight
355 1.0 that contains a list of vonMisesKappaConjugateRestraint on each
356 scale. If argument prior_rs is used, add them to that RestraintSet
359 if not (0 <= R <= c):
360 raise ValueError(
"parameters R and c should satisfy 0 <= R <= c")
364 prior_rs.add_restraint(
369 """scans the prot hierarchy and tries to find atom = (resno, name)
370 assumes that resno follows the same numbering as the sequence.
371 Stores already found atoms for increased speed.
373 if not hasattr(self,
'__memoized'):
374 self.__memoized = {prot: {}}
376 return self.__memoized[prot][atom]
381 residue_index=atom[0],
383 ).get_selected_particles()
385 print(
"found multiple atoms for atom %d %s!" % atom)
389 print(
"atom %d %s not found" % atom)
391 self.__memoized[prot][atom] = p0
396 Sets up a vonMises torsion angle restraint using the given kappa
397 particle as concentration parameter. Returns the restraint.
398 data is a list of observations.
404 Sets up a vonMises torsion angle restraint using the given kappa
405 particle as concentration parameter. Returns the restraint.
406 data is (mean, standard deviation).
408 raise NotImplementedError
412 Sets up a lognormal distance restraint using the given sigma and gamma.
413 Returns the restraint.
414 assumes atoms = (atom1, atom2)
415 where atom1 is (resno, atomname) and resno is the residue sequence
433 """Reads an ambiguous NOE restraint. contributions is a list of
434 (atom1, atom2) pairs, where atom1 is (resno, atomname). Sets up a
435 lognormal distance restraint using the given sigma and gamma.
436 Returns the restraint.
439 pairs = [(self.
find_atom(i, prot).get_particle(),
440 self.
find_atom(j, prot).get_particle())
for (i, j)
in
445 pairs, sigma, gamma, distance ** (-6))
449 self, prot, seqfile, tblfile, name=
'NOE', prior_rs=
None,
450 bounds_sigma=(1.0, 0.1, 100), bounds_gamma=(1.0, 0.1, 100),
451 verbose=
True, sequence_match=(1, 1)):
452 """read TBL file and store NOE restraints, using one sigma and one gamma
453 for the whole dataset. Creates the necessary uninformative priors.
454 - prot: protein hierarchy
455 - seqfile: a file with 3-letter sequence
456 - tblfile: a TBL file with the restraints
457 - name: an optional name for the restraintset
458 - prior_rs: when not None, add new sigma and gamma to this
459 RestraintSet instance.
460 - bounds_sigma or gamma: tuple of (initial value, lower, upper bound)
461 bounds can be -1 to set to default range [0,+inf]
462 - verbose: be verbose (default True)
463 - sequence_match : (noe_start, sequence_start)
464 Returns: data_rs, prior_rs, sigma, gamma
468 print(
"Prior for the NOE Scales")
476 first_residue_number=sequence_match[1])
477 tblr = IMP.isd.TBLReader.TBLReader(sequence,
478 sequence_match=sequence_match)
479 restraints = tblr.read_distances(tblfile,
'NOE')[
'NOE']
480 for i, restraint
in enumerate(restraints):
481 if verbose
and i % 100 == 0:
482 print(
"\r%d" % i, end=
' ')
487 if len(restraint[0]) > 1:
490 restraint[1], sigma, gamma)
493 restraint[1], sigma, gamma)
496 print(
"\r%d NOE restraints read" % i)
499 return rs, prior_rs, sigma, gamma
502 verbose=
True, sequence_match=(1, 1)):
503 """read TBL file and store NOE restraints, using the marginal of the
504 lognormal with one sigma and one gamma, for the whole dataset.
505 - prot: protein hierarchy
506 - seqfile: a file with 3-letter sequence
507 - tblfile: a TBL file with the restraints
508 - name: an optional name for the restraintset
509 - verbose: be verbose (default True)
510 - sequence_match : (noe_start, sequence_start)
517 first_residue_number=sequence_match[1])
518 tblr = IMP.isd.TBLReader.TBLReader(sequence,
519 sequence_match=sequence_match)
520 restraints = tblr.read_distances(tblfile,
'NOE')[
'NOE']
522 for i, restraint
in enumerate(restraints):
523 if verbose
and i % 100 == 0:
524 print(
"\r%d" % i, end=
' ')
530 pairs = [(self.
find_atom(i, prot).get_particle(),
531 self.
find_atom(j, prot).get_particle())
for (i, j)
in
534 ln.add_contribution(pairs, restraint[4])
537 print(
"\r%d NOE contributions added" % (len(restraints)))
542 """read TBL file and store lognormal restraints, using the marginal of the
543 lognormal with one sigma and gamma=1, for the whole dataset.
544 - prot: protein hierarchy
545 - seqfile: a file with 3-letter sequence
546 - tblfile: a TBL file with the restraints
547 - name: an optional name for the restraintset
548 - verbose: be verbose (default True)
555 tblr = IMP.isd.TBLReader.TBLReader(sequence)
556 restraints = tblr.read_distances(tblfile,
'HBond')[
'HBond']
558 for i, restraint
in enumerate(restraints):
559 if verbose
and i % 100 == 0:
560 print(
"\r%d" % i, end=
' ')
566 pairs = [(self.
find_atom(i, prot).get_particle(),
567 self.
find_atom(j, prot).get_particle())
for (i, j)
in
570 ln.add_contribution(pairs, restraint[4])
573 print(
"\r%d Hbond contributions added" % (len(restraints)))
577 sequence_match=(1, 1), name=
'TALOS', prior_rs=
None,
578 bounds_kappa=(1.0, 0.1, 10), verbose=
True, prior=
'jeffreys',
580 """read TALOS dihedral angle data, and create restraints for phi/psi
581 torsion angles, along with the prior for kappa, which is a scale for the
582 whole dataset, compare to 1/sigma**2 in the NOE case.
583 - prot: protein hierarchy
584 - seqfile: a file with 3-letter sequence
585 - talos_data: either a file (pred.tab), or a folder (pred/) in which all
586 files in the form res???.tab can be found. If possible,
587 try to provide the folder, as statistics are more
588 accurate in that case.
589 - fulldata : either True or False, whether the data is the full TALOS
590 output (predAll.tab or pred/ folder), or just the averages
592 - sequence_match : (talos_no, sequence_no) to adjust for different
594 - name: an optional name for the restraintset
595 - prior_rs: when not None, add new kappa(s) to this RestraintSet instance.
596 - bounds_kappa: tuple of (initial value, lower, upper bound)
597 bounds can be -1 to set to default range [0,+inf]
598 - verbose: be verbose (default True)
599 - prior: either 'jeffreys' or a tuple (R,c), which signifies to use the
600 conjugate prior of the von Mises restraint, with parameters R
601 and c. Good values are R=0 and c=10. Default: jeffreys prior.
602 - keep_all: in case of a folder for 'talos_data', whether to keep
603 candidates marked as 'outliers' by TALOS, or to include them.
604 Returns: data_rs, prior_rs, kappa
609 print(
"Prior for von Mises Kappa")
614 print(
"reading data")
617 first_residue_no=sequence_match[1])
620 sequence_match=sequence_match)
621 if os.path.isdir(talos_data):
623 for i, res
in enumerate(glob(os.path.join(talos_data,
'res???.tab'))):
624 if verbose
and i % 100:
625 print(
"\r%d" % i, end=
' ')
630 talosr.read(talos_data)
634 sequence_match=sequence_match)
635 talosr.read(talos_data)
637 data = talosr.get_data()
639 print(
"\rread dihedral data for %d residues" % len(data))
640 print(
"creating restraints")
642 for resno, datum
in data.items():
643 phidata = datum[
'phi']
644 psidata = datum[
'psi']
654 avgR.append(r.get_R0())
659 avgR.append(r.get_R0())
665 avgR.append(r.get_R0())
670 avgR.append(r.get_R0())
672 print(
"%s Restraints created. Average R0: %f" % \
673 (len(avgR), sum(avgR) / float(len(avgR))))
674 return rs, prior_rs, kappa
677 self, prot, profilefile, name=
'SAXS',
678 ff_type=IMP.saxs.HEAVY_ATOMS):
679 """read experimental SAXS profile and apply restraint the standard
681 Returns: a restraintset and the experimental profile
685 particles = IMP.atom.get_by_type(prot, IMP.atom.ATOM_TYPE)
687 rs.add_restraint(saxs_restraint)
688 return rs, saxs_profile
690 def _setup_md(self, prot, temperature=300.0, thermostat='berendsen',
691 coupling=500, md_restraints=
None, timestep=1.0, recenter=1000,
693 """setup molecular dynamics
694 - temperature: target temperature
695 - thermostat: one of 'NVE', rescale_velocities',
696 'berendsen', 'langevin'
697 - coupling: coupling constant (tau (fs) for berendsen,
698 gamma (/fs) for langevin)
699 - md_restraints: if not None, specify the terms of the energy to be used
700 during the md steps via a list of restraints.
701 - timestep: in femtoseconds.
702 - recenter: recenter the molecule every so many steps (Langevin only)
703 - momentum: remove angular momentum every so many steps (Berendsen only)
704 Returns: an instance of md and an instance of an OptimizerState (the
705 thermostat), or None if NVE.
709 md.assign_velocities(temperature)
710 md.set_time_step(timestep)
711 if thermostat ==
'NVE':
713 elif thermostat ==
'rescale_velocities':
716 md.add_optimizer_state(os)
717 elif thermostat ==
'berendsen':
720 md.add_optimizer_state(os)
723 md.add_optimizer_state(mom)
724 elif thermostat ==
'langevin':
727 md.add_optimizer_state(os)
732 raise NotImplementedError(thermostat)
735 md.set_restraints(md_restraints)
739 def _setup_normal_mover(self, particle, floatkey, stepsize):
740 """setup NormalMover to move particle's floatkey attribute
741 by a gaussian with standard deviation 'stepsize'
742 Returns: mover instance.
748 def _setup_md_mover(self, md, particles, temperature, n_md_steps=10):
749 """setup MDMover using md and particles.
751 - particles: particles to move, usually the leaves of the protein
753 - n_md_steps: number of md steps to perform on each move
754 Returns: mover instance.
758 return IMP.atom.MDMover(cont, md, temperature, n_md_steps)
760 def _setup_mc(self, mover, temperature=300.0, mc_restraints=None):
761 """setup monte carlo using a certain mover.
762 - mover: mover to use, NormalMover or MDMover usually.
763 - temperature: target temperature.
764 - mc_restraints: if not None, list of restraints for the metropolis
766 Returns: mc instance.
772 mc.set_kt(kB * temperature)
774 mc.set_return_best(
False)
777 mc.set_restraints(mc_restraints)
781 gamma=0.01, n_md_steps=10,
782 md_restraints=
None, mc_restraints=
None, timestep=1.0,
783 sd_threshold=0.0, sd_stepsize=0.01, sd_maxsteps=100):
784 """setup hybrid monte-carlo on protein. Uses basin hopping with steepest
785 descent minimization.
786 - prot: protein hierarchy.
787 - temperature: target temperature.
788 - gamma: coupling constant for langevin (/fs)
789 - n_md_steps: number of md steps per mc step
790 - md_restraints: if not None, specify the terms of the energy to be used
792 - mc_restraints: if not None, use these energy terms for the metropolis
794 - timestep: time step for md, in femtoseconds.
795 - sd_threshold: stop steepest descent after energy difference drops
797 - sd_stepsize: stepsize to use for the steepest descent, in angstrom.
798 - sd_maxsteps: maximum number of steps for steepest descent
799 Returns: hmc, mdmover, md and OptimizerState (thermostat)
801 raise NotImplementedError
802 md, os = self._setup_md(prot, temperature=temperature,
803 thermostat=
'langevin', coupling=gamma,
804 md_restraints=md_restraints, timestep=timestep)
805 particles = IMP.atom.get_by_type(prot, IMP.atom.ATOM_TYPE)
806 mdmover = self._setup_md_mover(md, particles, temperature, n_md_steps)
807 hmc = self._setup_mc(mdmover, temperature, mc_restraints)
809 sd.set_threshold(sd_threshold)
810 sd.set_step_size(sd_stepsize)
811 hmc.set_local_optimizer(sd)
812 hmc.set_local_steps(sd_maxsteps)
813 hmc.set_use_basin_hopping(
True)
814 return hmc, mdmover, md, os
817 n_md_steps=100, md_restraints=
None, mc_restraints=
None,
819 """setup hybrid monte-carlo on protein. Uses NVE MD and tries the full
820 - prot: protein hierarchy.
821 - temperature: target temperature.
822 - coupling: coupling constant (tau (fs) for berendsen,
823 gamma (/fs) for langevin)
824 - n_md_steps: number of md steps per mc step
825 - md_restraints: if not None, specify the terms of the energy to be used
827 - mc_restraints: if not None, use these energy terms for the metropolis
829 - timestep: time step for md, in femtoseconds.
830 - sd_threshold: stop steepest descent after energy difference drops
832 - sd_stepsize: stepsize to use for the steepest descent, in angstrom.
833 - sd_maxsteps: maximum number of steps for steepest descent
834 Returns: hmc, mdmover and md
836 prot = self._p[
'prot']
837 md, os = self._setup_md(prot, temperature=temperature,
838 thermostat=
'NVE', coupling=
None, timestep=1.0,
839 md_restraints=md_restraints)
840 particles = IMP.atom.get_by_type(prot, IMP.atom.ATOM_TYPE)
843 mdmover = PyMDMover(cont, md, n_md_steps)
849 mc.add_mover(mdmover)
851 mc.set_kt(kB * temperature)
853 mc.set_return_best(
False)
855 mc.set_move_probability(1.0)
856 return mc, mdmover, md
859 mc_restraints=
None, nm_stepsize=0.1):
860 """sets up monte carlo on nuisance, at a certain target temperature,
861 optionally using a certain set of restraints only.
862 - nuis: nuisance particle
863 - temperature: target temperature
864 - mc_restraints: optional set of restraints from which the energy should
865 be drawn instead of the energy of the complete system.
866 - floatkey: the floatkey to move.
867 - nm_stepsize: the stepsize of the normal mover
868 Returns: mc instance, nm instance.
870 nm = self._setup_normal_mover(nuis, nuis.get_nuisance_key(),
872 mc = self._setup_mc(nm, temperature, mc_restraints)
875 def _mc_and_update_nm(self, nsteps, mc, nm, stats_key,
876 adjust_stepsize=
True):
877 """run mc using a normal mover on a single particle,
878 update stepsize and store nsteps, acceptance and stepsize in the
879 statistics instance self.stat by using the key stats_key.
883 naccept = mc.get_number_of_forward_steps()
885 self.stat.increment_counter(stats_key, nsteps)
887 accept = float(naccept) / nsteps
888 self.stat.update(stats_key,
'acceptance', accept)
889 stepsize = nm.get_sigma()
890 self.stat.update(stats_key,
'stepsize', stepsize)
893 if 0.4 < accept < 0.6:
899 nm.set_sigma(stepsize * 2 * accept)
901 def _hmc_and_update_md(self, nsteps, hmc, mv, stats_key,
902 adjust_stepsize=
True):
903 """run hmc, update stepsize and print statistics. Updates number of MD
904 steps to reach a target acceptance between 0.4 and 0.6, sends
905 statistics to self.stat. MD steps are always at least 10 and at most 500.
908 naccept = hmc.get_number_of_forward_steps()
909 self.stat.increment_counter(stats_key, nsteps)
910 accept = float(naccept) / nsteps
911 self.stat.update(stats_key,
'acceptance', accept)
912 mdsteps = mv.get_number_of_steps()
913 self.stat.update(stats_key,
'n_md_steps', mdsteps)
915 if 0.4 < accept < 0.6:
917 mdsteps = int(mdsteps * 2 ** (accept - 0.5))
922 mv.set_nsteps(mdsteps)
925 """returns a string corresponding to the pdb structure of hierarchy
930 return output.getvalue()
932 def get_netcdf(self, prot):
933 raise NotImplementedError
936 temperature=300.0, prot_coordinates=
None):
937 """updates statistics for md simulation: target temp, kinetic energy,
938 kinetic temperature, writes coordinates and increments counter.
939 - md_key: stats md key
940 - nsteps: number of steps performed.
941 - md_instance: instance of the MolecularDynamics class.
942 - temperature: target temperature
943 - prot_coordinates: protein coordinates to be passed to the stats class,
946 self.stat.update(md_key,
'target_temp', temperature)
947 kinetic = md_instance.get_kinetic_energy()
948 self.stat.update(md_key,
'E_kinetic', kinetic)
949 self.stat.update(md_key,
'instant_temp',
950 md_instance.get_kinetic_temperature(kinetic))
951 self.stat.update_coordinates(md_key,
'protein', prot_coordinates)
952 self.stat.increment_counter(md_key, nsteps)
955 """shortcut for a frequent series of operations on MC simulations'
956 statistics. Creates an entry for acceptance, stepsize and one
957 coordinate set printed in the statistics file.
960 mc_key = stat.add_category(name=name)
962 stat.add_entry(mc_key, entry=Entry(
'temperature',
'%10G',
None))
963 stat.add_entry(mc_key, entry=Entry(
'acceptance',
'%10G',
None))
964 stat.add_entry(mc_key, entry=Entry(
'stepsize',
'%10G',
None))
966 stat.add_entry(mc_key, entry=Entry(coord,
'%10G',
None))
968 stat.add_entry(mc_key, name=
'counter')
972 """shortcut for a frequent series of operations on MD simulations'
973 statistics. Creates an entry for target temp, instantaneous temp,
974 kinetic energy, and one set of coordinates called 'protein' by
978 md_key = stat.add_category(name=name)
980 stat.add_entry(md_key, entry=Entry(
'target_temp',
'%10G',
None))
981 stat.add_entry(md_key, entry=Entry(
'instant_temp',
'%10G',
None))
982 stat.add_entry(md_key, entry=Entry(
'E_kinetic',
'%10G',
None))
984 stat.add_coordinates(md_key, coord)
986 stat.add_entry(md_key, name=
'counter')
990 """shortcut for a frequent series of operations on HMC simulations'
991 statistics. Adds acceptance, number of MD steps and a trajectory for
995 hmc_key = stat.add_category(name=name)
997 stat.add_entry(hmc_key, entry=Entry(
'temperature',
'%10G',
None))
998 stat.add_entry(hmc_key, entry=Entry(
'acceptance',
'%10G',
None))
999 stat.add_entry(hmc_key, entry=Entry(
'n_md_steps',
'%10G',
None))
1000 stat.add_entry(hmc_key, entry=Entry(
'E_kinetic',
'%10G',
None))
1002 stat.add_coordinates(hmc_key, coord)
1004 stat.add_entry(hmc_key, name=
'counter')
1008 """rescale the velocities of a bunch of particles having vx vy and vz
1014 p.set_value(k, p.get_value(k) * factor)
Applies a SingletonScore to each Singleton in a list.
def rescale_velocities
rescale the velocities of a bunch of particles having vx vy and vz floatkeys
CHARMMParameters * get_heavy_atom_CHARMM_parameters()
def init_model_ambiguous_NOE_restraint
Reads an ambiguous NOE restraint.
Maintains temperature during molecular dynamics.
Enforce CHARMM stereochemistry on the given Hierarchy.
Atoms get_phi_dihedral_atoms(Residue rd)
Return the atoms comprising the phi dihedral.
def do_md_protein_statistics
updates statistics for md simulation: target temp, kinetic energy, kinetic temperature, writes coordinates and increments counter.
def init_model_charmm_protein_and_ff
creates a CHARMM protein representation.
def find_atom
scans the prot hierarchy and tries to find atom = (resno, name) assumes that resno follows the same n...
def init_model_NOEs
read TBL file and store NOE restraints, using one sigma and one gamma for the whole dataset...
Various classes to hold sets of particles.
def init_simulation_setup_nuisance_mc
sets up monte carlo on nuisance, at a certain target temperature, optionally using a certain set of r...
A decorator for a particle which has bonds.
def init_simulation_setup_protein_hmc_nve
setup hybrid monte-carlo on protein.
Calculate score based on fit to SAXS profile.
A simple steepest descent optimizer.
void write_pdb(const Selection &mhd, TextOutput out, unsigned int model=1)
def init_model_TALOS
read TALOS dihedral angle data, and create restraints for phi/psi torsion angles, along with the prio...
Ambiguous NOE distance restraint between a number of pairs of particles.
def get_pdb
returns a string corresponding to the pdb structure of hierarchy prot.
def init_model_base
moves to wd and creates model
def init_model_NOEs_marginal
read TBL file and store NOE restraints, using the marginal of the lognormal with one sigma and one ga...
Return all close unordered pairs of particles taken from the SingletonContainer.
Classes to handle ISD statistics files.
def init_model_NOE_restraint
Sets up a lognormal distance restraint using the given sigma and gamma.
void read_pdb(TextInput input, int model, Hierarchy h)
CHARMM force field parameters.
Object used to hold a set of restraints.
Bond create_custom_bond(Bonded a, Bonded b, Float length, Float stiffness=-1)
Connect the two wrapped particles by a custom bond.
Class for storing model, its restraints, constraints, and particles.
Removes rigid translation and rotation from the particles.
def init_model_vonMises_restraint_full
Sets up a vonMises torsion angle restraint using the given kappa particle as concentration parameter...
Store a list of ParticleIndexPairs.
Conjugate prior for the concentration parameter of a von Mises distribution.
Maintains temperature during molecular dynamics by velocity scaling.
ParticleIndexPairs get_indexes(const ParticlePairsTemp &ps)
def init_model_jeffreys_kappa
given a list of scales, returns a RestraintSet('prior') with weight 1.0 that contains a list of vonMi...
Simple molecular dynamics optimizer.
Store a list of ParticleIndexes.
Atoms get_psi_dihedral_atoms(Residue rd)
Return the atoms comprising the psi dihedral.
def init_stats_add_md_category
shortcut for a frequent series of operations on MD simulations' statistics.
def init_model_HBonds_marginal
read TBL file and store lognormal restraints, using the marginal of the lognormal with one sigma and ...
Apply an NOE distance restraint between two particles.
Base class for all optimizers.
static Bonded setup_particle(Model *m, ParticleIndex pi)
Jeffreys prior for the concentration parameter of a von Mises distribution.
def init_stats_add_mc_category
shortcut for a frequent series of operations on MC simulations' statistics.
def init_simulation_setup_protein_hmc_hopper
setup hybrid monte-carlo on protein.
A decorator for a particle with x,y,z coordinates.
static Scale setup_particle(Model *m, ParticleIndex pi)
phi/psi dihedral restraint between four particles, using data from TALOS.
Score the bond based on a UnaryFunction,.
def read_sequence_file
read sequence of ONE chain, 1-letter or 3-letter, returns dict of no:3-letter code.
reads a TALOS file, or a TALOS folder, and stores the data
Modify a set of continuous variables using a normal distribution.
nonspecific methods used across all shared function objects.
void set_log_level(LogLevel l)
Set the current global log level.
A decorator for a residue.
Apply a lognormal distance restraint between two particles.
Basic functionality that is expected to be used by a wide variety of IMP users.
def init_model_standard_SAXS_restraint
read experimental SAXS profile and apply restraint the standard way (like foxs) Returns: a restraints...
def init_model_setup_scale
sets up a Scale particle to the initial default value.
static bool get_is_setup(const IMP::ParticleAdaptor &p)
def init_model_jeffreys
given a list of scales, returns a RestraintSet('prior') with weight 1.0 that contains a list of Jeffr...
Residue get_residue(Atom d, bool nothrow=false)
Return the Residue containing this atom.
A filter that excludes bonds, angles and dihedrals.
Class to handle individual particles of a Model object.
def init_model_vonMises_restraint_mean
Sets up a vonMises torsion angle restraint using the given kappa particle as concentration parameter...
def init_stats_add_hmc_category
shortcut for a frequent series of operations on HMC simulations' statistics.
def m
wrapper to call methods of m
Apply an NOE distance restraint between two particles.
Functionality for loading, creating, manipulating and scoring atomic structures.
CHARMMParameters * get_all_atom_CHARMM_parameters()
Hierarchies get_leaves(const Selection &h)
Classes to handle TALOS files or folders.
Select hierarchy particles identified by the biological name.
Applies a PairScore to each Pair in a list.
Support for small angle X-ray scattering (SAXS) data.
void set_check_level(CheckLevel tf)
Control runtime checks in the code.
def init_model_conjugate_kappa
given a list of scales, returns a RestraintSet('prior') with weight 1.0 that contains a list of vonMi...
Inferential scoring building on methods developed as part of the Inferential Structure Determination ...
Classes to handle TBL files.
Harmonic function (symmetric about the mean)
Maintains temperature during molecular dynamics.
A decorator for a particle with x,y,z coordinates and a radius.