- Contemporary messages sorted: [ by date ] [ by thread ] [ by subject ] [ by author ] [ by messages with attachments ]

From: wang <c00jsw00.nchc.org.tw>

Date: Thu, 03 May 2007 16:06:43 +0800

Dear all ,

I want to use distance constraint ( Ca atom of 121th residue and Ca of

251th residue ) .

Could you tell me how to do that ? I tried NOE module . But I still

didn't success .

thank you

Wang

David Mobley 提到:

*> It's probably also worth pointing out at this point that the Jarzynski
*

*> equality requires (on point 4) a large number (n) of simulations
*

*> because the most important values for the average are those at the
*

*> tails of the distribution. How large is hard to say -- but do some
*

*> testing at the beginning to see whether this is a price you can afford
*

*> to pay!
*

*>
*

*> (See the recent Jarzynski paper on convergence of exponentially
*

*> averaged work values). The Crooks approach is more efficient, but
*

*> requires doing both forward and reverse simulations for every
*

*> transformation. Not sure whether this is easy (possible?) to do with
*

*> the current code, as I'm not an SMD guy... You'd need fewer
*

*> realizations then.
*

*>
*

*> David
*

*>
*

*>
*

*> On 5/2/07, Gustavo Seabra <gustavo.seabra.gmail.com> wrote:
*

*>
*

*>> Hi James,
*

*>>
*

*>> What your code is doing, it seems to me, is averaging over all work
*

*>> values from *one* run, which has no meaning. To use the Jarzynski
*

*>> relationship, what you want is to, *at each point*, average a large
*

*>> number of work values. So, what you would do is:
*

*>>
*

*>> 1. Equilibrate the system, using a restraint to keep the coordinate of
*

*>> interest in the initial value. This will give you an equilibrated
*

*>> initial ensemble. (Now, if only we could define "equilibrate"... :-) )
*

*>>
*

*>> 2. From that equilibrated ensemble, pick a large number of structures
*

*>> (let's call it 'n'). These will be your initial structures for the next
*

*>> step, the pulling using the smd module.
*

*>>
*

*>> 3. (And this is the part you seem to be missing.) For each one of the
*

*>> 'n' structures you've picked, run an *independent* smd simulation,
*

*>> slowly pulling the coordinate from the initial to the final value.
*

*>> Notice that this will give you 'n' dist_vs_t files. Now the fun
*

*>> begins...
*

*>>
*

*>> 4. At every value of the reaction coordinate (1st column), average the
*

*>> work values from all the 'n' different simulations, using the Jarzynski
*

*>> formula. This will give you a (relative) free energy value at every
*

*>> value of the reaction coordinate, your "PMF" if you will.
*

*>>
*

*>> I hope that makes it a bit clearer for you.
*

*>>
*

*>> Good luck,
*

*>>
*

*>> Gustavo.
*

*>>
*

*>>
*

*>>
*

*>> James W wrote:
*

*>> > Dear Prof. Roitberg ,
*

*>> > I used the amber smd module to run a simulation . This is my output :
*

*>> > ____________________________________________________________________
*

*>> > x0(t) x force work
*

*>> > 57.00000 57.73576 -1471.51842 0.00000
*

*>> > 57.00093 57.34648 -691.10171 -1.00562
*

*>> > .....
*

*>> > ____________________________________________________________________
*

*>> > And this is my code for com[uting pmf :
*

*>> > ____________________________________________________________________
*

*>> > program pmf
*

*>> > real :: coor(100000),p(100000) ,work(100000)
*

*>> > beta = 1.04652438
*

*>> > open(unit=10,file='dist_vs_t')
*

*>> > do i=1,100000,1
*

*>> > read(10,*) a,b,c,d
*

*>> > coor(i) = a
*

*>> > work(i) = d
*

*>> > end do
*

*>> > ! PMF
*

*>> > do i=1,100000,1
*

*>> > pt = exp ( - work(i)/beta ) / i
*

*>> > p(i) = - beta * log (pt)
*

*>> > write(*,*) coor(i) , p(i)
*

*>> > end do
*

*>> > end program
*

*>> > ____________________________________________________________________
*

*>> > But the result is very curious , could you tell me that my code is
*

*>> right ?
*

*>> > thanks
*

*>> >
*

*>> > wang
*

*>> >
*

*>> >
*

*>> >
*

*>> >
*

*>> >
*

*>> >
*

*>> >
*

*>> >
*

*>> >
*

*>> >
*

*>> > On Wed, 02 May 2007 09:26:30 +0200, Adrian Roitberg wrote
*

*>> >
*

*>> >> James W wrote:
*

*>> >>
*

*>> >>> Dear all ,
*

*>> >>> #the reference : (http://en.wikipedia.org/wiki/Jarzynski_equality)
*

*>> >>> I confused the Jarzynski's equality . The equality is shown:
*

*>> >>> exp (-F/KT) = < exp (-W/KT) > .
*

*>> >>> &
*

*>> >>> F = < W > + ......
*

*>> >>> I used SMD module of CHARMM and obtained the reaction coordinates &
*

*>> >>> work ,like this :
*

*>> >>> ______________________________________________________
*

*>> >>> 57.00000 57.73576 -1471.51842 0.00000
*

*>> >>> 57.00093 57.34648 -691.10171 -1.00562
*

*>> >>> (RC) (WORK)
*

*>> >>> *RC : reaction coordinates
*

*>> >>> I wanted to obtain the " < W > " form my data . My method is :
*

*>> >>> 1. S = sum { exp (-W(i) /kT )} , i = 1 ....N
*

*>> >>> 2. p(i)= exp (-W(i) /kT ) / S
*

*>> >>> 3. <W(i)> = W(i)*p(i) / S
*

*>> >>>
*

*>> >>> Could you tell me that my method is right ?
*

*>> >>>
*

*>> >> James,
*

*>> >> I am not sure I fully understand your question, because you do not
*

*>> >> mention what the subindex i is in your formulas.
*

*>> >> Now, the formula from wikipedia is basically right, except for a
*

*>> couple
*

*>> >> of comments.
*

*>> >>
*

*>> >> First, it must be clear that F is really \Delta F. This is
*

*>> designed for
*

*>> >> free energy differences and NOT absolute free energies. Second, the
*

*>> >> exponential average of work values is over a set of runs, ALL
*

*>> >> starting from an equilibrated ensemble at a certain value of the
*

*>> >> coordinates (s) to be changed. The average then involves MANY runs.
*

*>> >>
*

*>> >> The second line F = < W > + ...... is a cumulant expansion of
*

*>> >> exponential average. I do not recommend using using this unless you
*

*>> >> really know what you are doing.
*

*>> >>
*

*>> >> Your formulas on how to implement this are not quite right.
*

*>> >> If you are thinking about the index i meaning time steps on one
*

*>> >> simulation, then it is not right.
*

*>> >> The subindex i should refer to the same distance and DIFFERENT
*

*>> simulations.
*

*>> >>
*

*>> >> Once that is done:
*

*>> >>
*

*>> >> 1. S = ( sum { exp (-W(i) /kT )}/N) , i = 1 ....N
*

*>> >> 2. \Delta F = - kT ln(S)
*

*>> >>
*

*>> >> I hope this makes things a bit clearer.
*

*>> >>
*

*>> >> Adrian
*

*>> >>
*

*>> -----------------------------------------------------------------------
*

*>> The AMBER Mail Reflector
*

*>> To post, send mail to amber.scripps.edu
*

*>> To unsubscribe, send "unsubscribe amber" to majordomo.scripps.edu
*

*>>
*

*> -----------------------------------------------------------------------
*

*> The AMBER Mail Reflector
*

*> To post, send mail to amber.scripps.edu
*

*> To unsubscribe, send "unsubscribe amber" to majordomo.scripps.edu
*

*>
*

-----------------------------------------------------------------------

The AMBER Mail Reflector

To post, send mail to amber.scripps.edu

To unsubscribe, send "unsubscribe amber" to majordomo.scripps.edu

Received on Sun May 06 2007 - 06:07:15 PDT

Date: Thu, 03 May 2007 16:06:43 +0800

Dear all ,

I want to use distance constraint ( Ca atom of 121th residue and Ca of

251th residue ) .

Could you tell me how to do that ? I tried NOE module . But I still

didn't success .

thank you

Wang

David Mobley 提到:

-----------------------------------------------------------------------

The AMBER Mail Reflector

To post, send mail to amber.scripps.edu

To unsubscribe, send "unsubscribe amber" to majordomo.scripps.edu

Received on Sun May 06 2007 - 06:07:15 PDT

Custom Search