[chimera-dev] implementing MMTK normal calculations in chimera
Conrad Huang
conrad at cgl.ucsf.edu
Tue Dec 8 17:59:00 PST 2009
Sorry to take so long to get back to you.
I think the problem is that MMTKinter creates MMTKChimeraModel instances
which derive from MMTK's ChemicalObjects.Molecule. Protein, on the
other hand, derives from ChemicalObjects.Complex. So there is no
compatible Protein object that you can pass to RigidMotionSubspace.
Looking at RigidMotionSubspace.__init__, it appears that it is expecting
a list of atom-containing objects as the second argument. So instead of
protein.residues(), I think you should be able to pass in a list of
MMTKChimeraGroup instances. That's pretty easy to construct from an
MMTKinter instance:
groups = []
for m in MMTKinter_instance.molecules:
groups.extend(m.groups)
The only trick is that protein.residues() only returns residues from
peptide chains. I'm not sure whether that's important for normal mode
calculation or was just used for convenience in the example. If it's
important to filter out non-amino acids, you will have to look at each
group. I _think_ that a group has attribute "peptide" if and only if it
is a standard amino acid.
Hope this helps.
Conrad
Victor Muñoz wrote:
> Dear developpers,
>
> I am trying to make an extension that allows normal mode calculations
> with MMTK into chimera. I have several problems beginning with the
> following.
>
> I have first used the MMTKinter.py and added a new function called
> NM_ffm that calculates normal modes with a full atom force field. This
> works perfectly (although I am not yet able to make an animation in
> chimera that follows the vectors calculated...that will be a future
> question for you). Then I added another function (NM_mcm) that allows
> me to do NM calculations with the mechanical constraint (see
> http://dirac.cnrs-orleans.fr/MMTK/using-mmtk/mmtk-example-scripts/normal-modes/normal-modes-with-mechanical-constraints).
> For this I need to define a subspace that requires a universe.protein
> object (the functions are given at the end of the mail). Unfortunately,
> I can not pass the molecular object from the chimera space as a
> universe.protein object. We thought in putting self.universe.protein
> call into the _makeUniverse function in MMTKinter.py but that did not
> work. We tried other combinations but we now face the fact that we are
> blocked. May be we need another function in order to define subspaces
> but yet we still don't know how to pass chimera infos to MMTK.
>
> Could you give me a hand so that I can define correctly subspaces on
> chimera molecules?
>
> All the best and thank in advance for any help,
>
> Victor
>
> *************
> Computational Biotechnological Chemistry @transmet
>
>
>
>
>
> ####Function for normal modes calculations in full atom ff
> def NM_ffm(self, nsteps, stepsize=0.02,
> interval=None, action=None, **kw):
> from chimera import replyobj
> timestamp("_minimize")
> from MMTK import Units
> from MMTK.ForceFields.Amber import AmberData
> from MMTK.Minimization import SteepestDescentMinimizer
> from MMTK.Trajectory import LogOutput
> # from JD to make NM calculations
> from MMTK.NormalModes import NormalModes
> from chimera import runCommand as run
> from chimera import selection as sele
> from chimera import Point, Vector
> import os, sys
> # JD'end
> import sys
> if not interval:
> actions = []
> else:
> actions = [ LogOutput(sys.stdout, ["energy"],
> interval, None, interval) ]
> kw["step_size"] = stepsize * Units.Ang
> minimizer = SteepestDescentMinimizer(self.universe,
> actions=actions, **kw)
> if action is None or not interval:
> interval = None
> msg = "Initial energy: %f" % self.universe.energy()
> replyobj.status(msg)
> replyobj.info <http://replyobj.info>(msg)
> saveNormalizeName = AmberData._normalizeName
> AmberData._normalizeName = simpleNormalizeName
> remaining = nsteps
> while remaining > 0:
> timestamp(" minimize interval")
> if interval is None:
> realSteps = remaining
> else:
> realSteps = min(remaining, interval)
> minimizer(steps=realSteps)
> remaining -= realSteps
> if action is not None:
> action(self)
> timestamp(" finished %d steps" % realSteps)
> msg = "Finished %d of %d minimization steps" % (
> nsteps - remaining, nsteps)
> replyobj.status(msg)
> replyobj.info <http://replyobj.info>(msg)
> replyobj.info <http://replyobj.info>("\n")
> #Calculate normal modes JD's trucho
> def despl(vec):
> run("sel #0")
> atoms=sele.currentAtoms()
> for a in atoms:
> a.setCoord(a.xformCoord()-10*(vec))
> modes = NormalModes(self.universe)
> for mode in modes:
> print mode
> print mode.array
> L1=[]
> L2=[]
> L3=[]
> for mode in modes:
> y=6
> u=len(mode.array)
> while y<u:
> L1.append(mode.array[y][0])
> L2.append(mode.array[y][1])
> L3.append(mode.array[y][2])
> y+=1
> y=0
> run("sel:all")
> atoms=sele.currentAtoms()
> run("~sel")
> u=len(L1)
> a=0
> o=len(atoms)
> while y<u:
> vec=Vector(float(L1[y]),float(L2[y]),float(L3[y]))
> atoms[a].setCoord(atoms[a].coord()-(vec))
> a+=1
> if a==o:
> a=0
> raw_input()
> else:
> pass
> y+=1
> #JD's end
> return modes
> AmberData._normalizeName = saveNormalizeName
> timestamp("end _minimize")
>
> ####Function for mechanical constraint
> def NM_mcm(self, nsteps, stepsize=0.02,
> interval=None, action=None, **kw):
> from MMTK import *
> from MMTK.Proteins import Protein
> from MMTK.ForceFields import Amber94ForceField
> from MMTK.NormalModes import NormalModes, SubspaceNormalModes
> from MMTK.Subspace import RigidMotionSubspace
> from MMTK.Minimization import ConjugateGradientMinimizer
> from MMTK.Trajectory import StandardLogOutput
> import numpy
> # Victor Modifications
> from chimera import replyobj
> timestamp("_minimize")
> from MMTK import Units
> from MMTK.ForceFields.Amber import AmberData
> from MMTK.Minimization import SteepestDescentMinimizer
> from MMTK.Trajectory import LogOutput
> from MMTK.NormalModes import NormalModes
> from chimera import runCommand as run
> from chimera import selection as sele
> from chimera import Point, Vector
> import os, sys
> #/end
> if not interval:
> actions = []
> else:
> actions = [ LogOutput(sys.stdout, ["energy"],
> interval, None, interval) ]
> kw["step_size"] = stepsize * Units.Ang
> minimizer = SteepestDescentMinimizer(self.universe,
> actions=actions, **kw)
> if action is None or not interval:
> interval = None
> msg = "Initial energy: %f" % self.universe.energy()
> replyobj.status(msg)
> replyobj.info <http://replyobj.info>(msg)
> saveNormalizeName = AmberData._normalizeName
> AmberData._normalizeName = simpleNormalizeName
> remaining = nsteps
> while remaining > 0:
> timestamp(" minimize interval")
> if interval is None:
> realSteps = remaining
> else:
> realSteps = min(remaining, interval)
> minimizer(steps=realSteps)
> remaining -= realSteps
> if action is not None:
> action(self)
> timestamp(" finished %d steps" % realSteps)
> msg = "Finished %d of %d minimization steps" % (
> nsteps - remaining, nsteps)
> replyobj.status(msg)
> replyobj.info <http://replyobj.info>(msg)
> replyobj.info <http://replyobj.info>("\n")
> # Construct system
> #self.universe = InfiniteUniverse(Amber94ForceField())
> #self.universe.protein = self.universe
> # Set up the subspace: rigid-body translation and rotation for
> each residue
> self.subspace = RigidMotionSubspace(self.universe,
> self.universe.protein.residues())
>
>
>
>
> ------------------------------------------------------------------------
>
> _______________________________________________
> Chimera-dev mailing list
> Chimera-dev at cgl.ucsf.edu
> http://www.cgl.ucsf.edu/mailman/listinfo/chimera-dev
More information about the Chimera-dev
mailing list