# Re: FEP question

From: luca (bellucci14_at_unisi.it)
Date: Fri Sep 07 2007 - 12:43:24 CDT

Hi Jerome
I am in agreement with you.
With inifinite masses I meant very big masses like 1E40.
In this way ( with initial velocity equal zero!!!)
you can make same behavior of fixed atoms.....

Luca

> Hi Luca,
> Fixed atoms do constitute a bias. I'm afraid the argument about
> infinite masses, though formally interesting, does not change that
> state of affairs.
> One can adopt either point of view:
> 1) Constraints are applied by adding constraint forces (deriving from
> a constrained Hamiltonian). Hence a bias in configurational sampling.
> 2) The Hamiltonian is unconstrained, but constraints are applied by
> setting infinite masses on a set of atoms. As a consequence, that part
> of the system has an infinite relaxation time. Obviously, relaxation
> does not occur during the simulation, which is therefore not sampling
> the equilibrium ensemble, but another ensemble dictated by the initial
> conditions, aka the constrained ensemble.
>
> So the "no force, no bias" idea only goes so far...
> Jerome
>
> On 9/7/07, luca <bellucci14_at_unisi.it> wrote:
> > Hi all,
> > I think you are right....but I have one observation.
> > Can you test use of fixed atom only in initial step?
> > In this way you do not have any additional force therefore no bias.
> > If you set fixed atom you have put infinite mass to some atoms but you do not
> > change potential form of hamiltonian value.
> > In this way you Infuence only the sampling at initial step.
> > However you should minimize the number of restraint for a good sampling.
> >
> > Best,
> > Luca
> >
> >
> > > Matthew,
> > > I think your ideas are essentially the right things to do.
> > > In principle, you can use a custom Tcl script that enforces soft
> > > restraints to, as you said, keep the two rings roughly in the same
> > > place.
> > >
> > > When dealing with isolated molecules, this is not a problem at all,
> > > since at the end-points (lambda=0 or 1), one of the molecules is
> > > virtually non-existent, so tethering it to the other does not bias the
> > > energies. Only points at 0 < lambda < 1 are affected, which just means
> > > that the free energy is computed along a different pathway - a
> > > better-behaved pathway, actually.
> > > Now, the issue with your DNA residues is that even at the end-points,
> > > the ghost residue still has bonded terms, so tethering it may induce a
> > > bias.
> > >
> > > Still, I see two possible courses of action:
> > > 1) add such restraints coupling the two residues, but in a way that
> > > does not add a significant conformational bias - this may be possible
> > > if "native" conformational fluctuations are small.
> > > 2) as you also suggest, add torsional restraints (same kind of custom
> > > Tcl script) to the ghost residue until its nonbonded interactions are
> > > strong enough. The issue there is that in theory, it would require an
> > > FEP step to determine the FE difference between the restrained residue
> > > and the unrestrained one... this may be done by growing the restraints
> > > using the "conformational free energy" module of NAMD, but it is still
> > > adding a nontrivial step to a nontrivial calculation.
> > >
> > > You see, there is no cheap way out, but I hope it helps nonetheless.
> > > Best,
> > > Jerome
> > >
> > >
> > > On 9/6/07, Matthew WIlce <matthew.wilce_at_med.monash.edu.au> wrote:
> > > > Dear NAMDers,
> > > > I have a couple of questions regarding FEP.
> > > >
> > > > I am using FEP and mutating a Cytosine into a Thymine, ADE and GUA. I
> > > > have set up the dual topologies so that it is the purine or
> > > > pyrimidine that is mutating while the rest of residue in not
> > > > changing. I have noticed that at the beginning of the FEP run the
> > > > purine or pyrimidine that is starting to grow moves around too much
> > > > and sometimes the ring flips ~180 degrees. When the ring flips 180
> > > > degrees it prevents a sensible outcome from the FEP. This also occurs
> > > > for the Cytosine ring towards the end of the FEP. Once lambda has
> > > > become greater than about 0.01 and the growing ring starts to
> > > > interact with the rest of the molecules it behaves much better.
> > > >
> > > > For a purine to purine mutation eg ADE to GUA I could of course
> > > > modify the dual topology so that the all of the common atoms are
> > > > maintained throughout the FEP, and the same goes for the pyrimidines.
> > > > But I cannot see how a work around when the rings are quite different
> > > > eg. pyrimidine vs purine.
> > > >
> > > > Is there a way to sensibly restrain/constrain the movement of the
> > > > growing component of the FEP so that it cannot 'float' around until
> > > > lambda1 passes say 0.01, and also apply similar restraints to the
> > > > shrinking component when lambda 1 is greater then ~0.99.
> > > >
> > > > Has anyone done this? What impact could this have on the free energies?
> > > >
> > > > I have thought of writing a tcl script that constrains the atoms at
> > > > the beginning and end of the FEP depending on lambda ?? but I am not
> > > > sure how. I want the growing and shrinking rings of the purine or
> > > > pyrimidine to stay roughly in the same place relative to each other
> > > > until their interactions with the rest of the molecule may cause them
> > > > to move if necessary and not have them moving because they are not
> > > > interacting.
> > > >
> > > > I have tried using the namd parameters as suggested in the most
> > > > recent FEP tutorial as well as varying them but I have not been able
> > > > to resolve this. Even using a very low number of steps in each FEP
> > > > window (as low as 1000) does not stop the ring from sometimes flipping.
> > > >
> > > > Thanks in advance,
> > > > Matthew
> > > >
> > > >
> > >
> > >
> >
> >
> >
>

This archive was generated by hypermail 2.1.6 : Wed Feb 29 2012 - 15:45:13 CST