Re: Doubts and problem in ABF simulations

From: Giacomo Fiorin (gfiorin_at_seas.upenn.edu)
Date: Wed Sep 02 2009 - 13:21:09 CDT

Hi Paulo, it definitely looks like, in some cases, we should allow the
colvars module to handle the dummy atom with ABF as well.

However, in your application: are you sure that the distance from the
dummy atom is more appropriate? Unless you fix the protein at an
absolute position in space, the colvar you just defined will be the
distance from the origin (or whatever point in space you set after
dummyAtom), and NOT from the binding site of the protein.

Also, take into account that when using SMD, you steer the ligand
along a specific trajectory. When using ABF over a distance variable,
you also sample movements of the ligand orthogonal to that trajectory.
 So the phase spaces covered with the two simulations will be
different.

Giacomo

---- ----
 Giacomo Fiorin
   ICMS - Institute for Computational Molecular Science
     Temple University
     1900 N 12 th Street, Philadelphia, PA 19122
 work phone: (+1)-215-204-4216
 mobile: (+1)-267-324-7676
 mail: giacomo.fiorin_at_gmail.com
---- ----

On Wed, Sep 2, 2009 at 11:55 AM, Paulo Cesar Telles de
Souza<paulocts_at_gmail.com> wrote:
> Dear Jerome
>
> Thanks a lot for your help! I didn't know that the NAMD 2.7b1 altered ABF
> method (and others) in such a way. I have decided to use this new version of
> NAMD with colvars module. I don't think that it is necessary employing
> multidimensional ABF for my problem. I have obtained PMF using a SMDs and
> Jarzynski equality and now I would like to compare the result with one
> obtained using another free energy calculation method. To specify the same
> reaction coordinate that I have used in SMD, I have tried to use the
> distance component, that is, I chose ligand atoms for group 1 and dummyAtom
> option for group 2 (setting the appropriate coordinate to reproduce the same
> reaction coordinate), but, still, it didn't work. I pasted the error message
> below. Is it necessary to use another option together with dummyAtom?
>
> Bests,
> Paulo
>
> colvars:
> ----------------------------------------------------------------------
> colvars: Initializing the collective variables module, version 29-09-2008.
> colvars: # colvarsTrajFrequency = 500
> colvars: # colvarsRestartFrequency = 500
> colvars: # trajAppend = off [default]
> colvars: The restart output state file will be "protabf-a.colvars.state".
> colvars: The trajectory file will be "protabf-a.colvars.traj".
> colvars: The final output state file will be "protabf-a.colvars.state".
> colvars: # analysis = off [default]
> colvars:
> ----------------------------------------------------------------------
> colvars: Initializing a new collective variable.
> colvars: # name = Distance
> colvars: Initializing a new "distance" component.
> colvars: # componentCoeff = 1 [default]
> colvars: # componentExp = 1 [default]
> colvars: # oneSiteSystemForce = off [default]
> colvars: Initializing atom group "group1".
> colvars: Atom group "group1" defined, 35 atoms, total mass = 650.977.
> colvars: Initializing atom group "group2".
> colvars: Atom group "group2" defined, 0 atoms, total mass = 1.
> colvars: All components initialized.
> colvars: # width = 0.1
> colvars: # lowerBoundary = 16
> colvars: Lower boundary defined.
> colvars: # upperBoundary = 31
> colvars: Upper boundary defined.
> colvars: # lowerWallConstant = 10
> colvars: Applying a harmonic lower wall at 16, with coefficient 10.
> colvars: # upperWallConstant = 10
> colvars: Applying a harmonic upper wall at 31, with coefficient 10.
> colvars: # extendedLagrangian = off [default]
> colvars: # outputValue = on [default]
> colvars: # outputVelocity = off [default]
> colvars: # outputSystemForce = off [default]
> colvars: # outputAppliedForce = off [default]
> colvars: ------------------------------
> ----------------------------------------
> colvars: Collective variables initialized, 1 in total.
> colvars:
> ----------------------------------------------------------------------
> colvars: Initializing a new "abf" instance.
> colvars: # name = "abf1" [default]
> colvars: # colvars = { Distance }
> colvars: # applybias = on [default]
> colvars: # hidejacobian = off [default]
> colvars: Jacobian (geometric) forces will be included in reported free
> energy gradients.
> colvars: # fullsamples = 500
> colvars: # inputprefix = [default]
> colvars: # outputfreq = 500 [default]
> colvars: Finished ABF setup.
> colvars:
> ----------------------------------------------------------------------
> colvars: Collective variables biases initialized, 1 in total.
> colvars:
> ----------------------------------------------------------------------
> colvars: Collective variables module initialized.
> colvars:
> ----------------------------------------------------------------------
> Info: Startup phase 7 took 0.012846 s, 7.72452 MB of memory in use
> Info: Startup phase 8 took 0.000380039 s, 12.3929 MB of memory in use
> Info: Finished startup at 0.541487 s, 12.3929 MB of memory in use
>
> TCL: Running for 4000000 steps
> ETITLE: TS BOND ANGLE DIHED
> IMPRP ELECT VDW BOUNDARY
> MISC KINETIC TOTAL TEMP
> POTENTIAL TOTAL3 TEMPAVG PRESSURE
> GPRESSURE VOLUME PRESSAVG GPRESSAVG
>
> ENERGY: 0 849.8375 2408.4881 1348.5871
> 135.3159 -5237.1600 -1091.2560 0.0000
> 0.0000 3104.3467 1518.1595 292.6250 -1586.1873
> 1557.9272 292.6250 112.1152 111.0410
> 590148.0000 112.1152 111.0410
>
> colvars: Error: system forces are not available from a dummy atom group.
> colvars: If this error message is unclear, try recompile with
> -DCOLVARS_DEBUG.
> FATAL ERROR: Error in the collective variables module: exiting.
>
> ------------- Processor 0 Exiting: Called CmiAbort ------------
> Reason: FATAL ERROR: Error in the collective variables module: exiting.
>
>
> [0] Stack Traceback:
> [0] CmiAbort+0x7f [0xa4526d]
> [1] _Z8NAMD_diePKc+0x62 [0x500082]
> [2] _ZN16colvarproxy_namd11fatal_errorERKSs+0x64 [0x844b68]
> [3] _ZN12colvarmodule11fatal_errorERKSs+0x18 [0x79ae22]
> [4] _ZNK12colvarmodule10atom_group12system_forceEv+0x53 [0x816f75]
> [5] _ZN6colvar8distance19calc_force_invgradsEv+0xd7 [0x829145]
> [6] _ZN6colvar4calcEv+0xc03 [0x7cca1f]
> [7] _ZN12colvarmodule4calcEv+0x6c [0x7969bc]
> [8] _ZN16colvarproxy_namd9calculateEv+0x497 [0x8415a1]
> [9] _ZThn16_N16colvarproxy_namd9calculateEv+0xa [0x841108]
> [10]
> _ZN12GlobalMaster11processDataEPiS0_P6VectorS2_S2_S0_S0_S2_S0_S0_S2_+0x71
> [0x77b165]
> [11] _ZN18GlobalMasterServer11callClientsEv+0x469 [0x77f191]
> [12] _ZN18GlobalMasterServer8recvDataEP20ComputeGlobalDataMsg+0x6f7
> [0x77e7d7]
> [13] _ZN10ComputeMgr21recvComputeGlobalDataEP20ComputeGlobalDataMsg+0x12
> [0x569388]
> [14]
> _ZN18CkIndex_ComputeMgr48_call_recvComputeGlobalData_ComputeGlobalDataMsgEPvP10ComputeMgr+0xf
> [0x569373]
> [15] CkDeliverMessageFree+0x21 [0x9c2d71]
> [16] _Z15_processHandlerPvP11CkCoreState+0x509 [0x9c2365]
> [17] CsdScheduleForever+0xa5 [0xa4bcf5]
> [18] CsdScheduler+0x1c [0xa4b8f6]
> [19] _ZN7BackEnd7suspendEv+0xb [0x508add]
> [20] _ZN9ScriptTcl7Tcl_runEPvP10Tcl_InterpiPPc+0x140 [0x8fcde4]
> [21] TclInvokeStringCommand+0x91 [0xa6e208]
> [22] namd2 [0xaa4058]
> [23] Tcl_EvalEx+0x176 [0xaa469b]
> [24] Tcl_EvalFile+0x134 [0xa9c0a4]
> [25] _ZN9ScriptTcl3runEPc+0x14 [0x8fc4e2]
> [26] _Z18after_backend_initiPPc+0x22b [0x50481b]
> [27] main+0x3a [0x5045ba]
> [28] __libc_start_main+0xe6 [0x7ff4b89e8466]
> [29] _ZNSt8ios_base4InitD1Ev+0x52 [0x4ff9ea]
> Charm++ fatal error:
> FATAL ERROR: Error in the collective variables module: exiting.
>
>
> Aborted
> 2009/8/31 Jerome Henin <jhenin_at_cmm.chem.upenn.edu>
>>
>> Dear Paulo,
>>
>> About your doubt: it is not easy to describe ligand-receptor binding
>> with a single parameter. Sop far, I am not aware of anyone
>> successfully applying ABF to a protein-ligand binding problem. You may
>> have better luck with a multidimensional calculation, using either ABF
>> or metadynamics.
>>
>> About your problem: indeed, there seems to be an error in the
>> documentation for the abscissa parameter. As mentioned prominently on
>> the webpage you link to, the NAMD 2.6 ABF code is deprecated, and only
>> the colvars module of NAMD 2.7b1 is actively maintained. If you cannot
>> use NAMD 2.7b1, one option is to modify the abscissa.tcl script to
>> accept an atom group for abf3 rather than a single atom. This should
>> be relatively straightforward. However, it will not solve the problem
>> of choosing the right coordinate(s) for the job.
>>
>> Best,
>> Jerome
>>
>>
>> On Mon, Aug 31, 2009 at 12:00 PM, Paulo Cesar Telles de
>> Souza<paulocts_at_gmail.com> wrote:
>> > Dear all,
>> >
>> > I have some doubts and one problem on running ABF simulations in NAMD.
>> >
>> > Doubts: I'm not sure what order parameter is better for my system,
>> > mainly
>> > because I have tested to use the z-coord order parameter and it didn't
>> > work.
>> > I selected some atoms in the binding site (abf1) and the ligant atoms
>> > (abf2). The z axis of the box is parallel to the direction from the
>> > ligand
>> > in the receptor to the solvent bulk and the distance vector between abf1
>> > and
>> > abf2 is almost parallel to the z-direction. For my suprise, when I
>> > performed
>> > the ABF simulation, the ligant "walk'' in the protein surface. Did I
>> > choose
>> > the wrong order parameter or the problem is due to my vector abf1-abf2 ?
>> > Does the vector abf1-abf2 change during ABF simulations ?
>> >
>> > Problem: I have tried to use the abscissa order parameter to obtain a
>> > PMF
>> > along a dissociation pathway through which a ligand unbinds from its
>> > receptor (the reaction coordinate is ligant-protein to ligant-solvent).
>> > I
>> > used a group of atoms in abf3 but received the following message error.
>> > This
>> > is strange because the manual and this repository
>> > (http://www.edam.uhp-nancy.fr/ABF/repository.html) specify that I can
>> > choose
>> > a list of atoms. If I try only with one atom, it works. How can I use a
>> > list
>> > in abf3 ( I want that a entire system to move along the vector defined
>> > by
>> > abf1 and abf2) ?
>> >
>> > Thanks on advance
>> >
>> > Paulo
>> >
>> > TCL: ABF> ---------------------------------------------
>> > TCL: ABF> Adaptive Biasing Force protocol version 1.7
>> > TCL: ABF> ---------------------------------------------
>> > TCL: ABF>
>> > TCL: ABF> Using coordinate type : abscissa
>> > TCL: ABF> Position of an atom along the axis joining two groups
>> > TCL: ABF> dxi : 0.1
>> > TCL: ABF> outFile : abf.dat [default]
>> > TCL: ABF> fullSamples : 500
>> > TCL: ABF> writeFxiFreq : 0 [default]
>> > TCL: ABF> writeXiFreq : 100
>> > TCL: ABF> abf2 : 4240 4241 4242 4243 4244 4245
>> > TCL: ABF> fMax : 60.0 [default]
>> > TCL: ABF> outputFreq : 5000 [default]
>> > TCL: ABF> inFiles : [default]
>> > TCL: ABF> historyFile : none [default]
>> > TCL: ABF> df : 1.0 [default]
>> > TCL: ABF> moveBoundary : 0 [default]
>> > TCL: ABF> distFile : none [default]
>> > TCL: ABF> temp : 300.0 [default]
>> > TCL: ABF> forceConst : 10.0
>> > TCL: ABF> dSmooth : 0.1
>> > TCL: ABF> xiMin : -1.0
>> > TCL: ABF> xiMax : 15.0
>> > TCL: ABF> abf1 : 2724 2725 2726 2727 2728 2729 2730 2731
>> > 2732
>> > 2733 2734 2735 2736 2737 2738 2739 2806 2807 2808 2809 2810 2811 2812
>> > 2813
>> > 2814 2815 2816 2817 2818 2819 2820 2821 2822 2823 2824 2825 2826 2827
>> > 2828
>> > 2829 2830 2831 3040 3041 3042 3043 3044 3045 3046 3047 3048 3049 3050
>> > 3051
>> > 3052 3053 3054 3055 3056 3103 3104 3105 3106 3107 3108 3109 3110 3111
>> > 3112
>> > 3113 3114 3115 3116 3117 3118 3119 3120 3121 3418 3419 3420 3421 3422
>> > 3423
>> > 3424 3425 3426 3427 3428 3429 3430 3431 3432 3433 3434 3435 3436 3437
>> > 3438
>> > 3439 3440 3441 3497 3498 3499 3500 3501 3502 3503 3504 3505 3506 3507
>> > 3508
>> > 3509 3510 3511 3512 3513 3514 3515
>> > TCL: ABF> Accumulating force data in 160 bins
>> > TCL: ABF> Moving abf.dat to abf.dat.BAK
>> > TCL: ABF> Using the abscissa of atom 4240 4243
>> > TCL: expected integer but got "4240 4243"
>> > FATAL ERROR: expected integer but got "4240 4243"
>> > while executing
>> > "addatom $abf3"
>> > (in namespace eval "::ABF::ABFcoord" script line 9)
>> > invoked from within
>> > "namespace eval ABFcoord {
>> >
>> > set abf1 $::ABF::abf1
>> > set abf2 $::ABF::abf2
>> > set abf3 $::ABF::abf3
>> >
>> > print "ABF> Using the abscissa of atom $abf3"
>> >
>> > adda..."
>> > (procedure "ABFstartup" line 2)
>> > invoked from within
>> > "ABFstartup"
>> > (in namespace eval "::ABF" script line 4)
>> > invoked from within
>> > "namespace eval ::ABF {
>> >
>> > # coordinate-specific startup procedure
>> > ABFstartup
>> >
>> > if { [info exists restraintList] } {
>> >
>> > array set rArray $restraintList
>> > ..."
>> > (file "/home/paulo/Desktop/Tutorial-ABF/abf-1.8/abf_script.tcl" line
>> > 766)
>> > ------------- Processor 0 Exiting: Called CmiAbort ------------
>> > Reason: FATAL ERROR: expected integer but got "4240 4243"
>> > while executing
>> > "addatom $abf3"
>> > (in namespace eval "::ABF::ABFcoord" script line 9)
>> > invoked from within
>> > "namespace eval ABFcoord {
>> >
>> > set abf1 $::ABF::abf1
>> > set abf2 $::ABF::abf2
>> > set abf3 $::ABF::abf3
>> >
>> > print "ABF> Using the abscissa of atom $abf3"
>> >
>> > adda..."
>> > (procedure "ABFstartup" line 2)
>> > invoked from within
>> > "ABFstartup"
>> > (in namespace eval "::ABF" script line 4)
>> > invoked from within
>> > "namespace eval ::ABF {
>> >
>> > # coordinate-specific startup procedure
>> > ABFstartup
>> >
>> > if { [info exists restraintList] } {
>> >
>> > array set rArray $restraintList
>> > ..."
>> > (file "/home/paulo/Desktop/Tutorial-ABF/abf-1.8/abf_script.tcl" line
>> > 766)
>> >
>> > [0] Stack Traceback:
>> > [0] CmiAbort+0x7f [0xa4526d]
>> > [1] _Z8NAMD_diePKc+0x62 [0x500082]
>> > [2] _ZN15GlobalMasterTcl10initializeEv+0x130 [0x789f22]
>> > [3] _ZN15GlobalMasterTclC9Ev+0x66 [0x786426]
>> > [4] _ZN15GlobalMasterTclC1Ev+0x6 [0x78646c]
>> > [5] _ZN10ComputeMgr14createComputesEP10ComputeMap+0x103 [0x56cb59]
>> > [6] _ZN4Node7startupEv+0x287 [0x89f303]
>> > [7] _ZN12CkIndex_Node18_call_startup_voidEPvP4Node+0x12 [0x89f078]
>> > [8] CkDeliverMessageFree+0x21 [0x9c2d71]
>> > [9] _Z15_processHandlerPvP11CkCoreState+0x509 [0x9c2365]
>> > [10] CsdScheduleForever+0xa5 [0xa4bcf5]
>> > [11] CsdScheduler+0x1c [0xa4b8f6]
>> > [12] _ZN7BackEnd7suspendEv+0xb [0x508add]
>> > [13] _ZN9ScriptTcl9initcheckEv+0x80 [0x90041e]
>> > [14] _ZN9ScriptTcl7Tcl_runEPvP10Tcl_InterpiPPc+0x20 [0x8fccc4]
>> > [15] TclInvokeStringCommand+0x91 [0xa6e208]
>> > [16] /home/paulo/programas/NAMD/namd2 [0xaa4058]
>> > [17] Tcl_EvalEx+0x176 [0xaa469b]
>> > [18] Tcl_EvalFile+0x134 [0xa9c0a4]
>> > [19] _ZN9ScriptTcl3runEPc+0x14 [0x8fc4e2]
>> > [20] _Z18after_backend_initiPPc+0x22b [0x50481b]
>> > [21] main+0x3a [0x5045ba]
>> > [22] __libc_start_main+0xe6 [0x7ff1a10d4466]
>> > [23] _ZNSt8ios_base4InitD1Ev+0x52 [0x4ff9ea]
>> > Fatal error on PE 0> FATAL ERROR: expected integer but got "4240 4243"
>> > while executing
>> > "addatom $abf3"
>> > (in namespace eval "::ABF::ABFcoord" script line 9)
>> > invoked from within
>> > "namespace eval ABFcoord {
>> >
>> > set abf1 $::ABF::abf1
>> > set abf2 $::ABF::abf2
>> > set abf3 $::ABF::abf3
>> >
>> > print "ABF> Using the abscissa of atom $abf3"
>> >
>> > adda..."
>> > (procedure "ABFstartup" line 2)
>> > invoked from within
>> > "ABFstartup"
>> > (in namespace eval "::ABF" script line 4)
>> > invoked from within
>> > "namespace eval ::ABF {
>> >
>> > # coordinate-specific startup procedure
>> > ABFstartup
>> >
>> > if { [info exists restraintList] } {
>> >
>> > array set rArray $restraintList
>> > ..."
>> > (file "/home/paulo/Desktop/Tutorial-ABF/abf-1.8/abf_script.tcl" line
>> > 766)
>
>

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