[gmx-users] Free energy perturbation parameters sc-alpha sc-power sc-sigma

bipin singh bipinelmat at gmail.com
Thu Jan 23 10:56:41 CET 2014


Ok, now I understand the transformation scheme. I think for doing the
mentioned transformation you have to do the following things:

1) In the first set of transformation, transform the HG1 (1,2 and 3) and
HG2 (1,2 and 3) to dummy atoms (you have written "I have to turn on the HG"
can you explain why?) by first relaxing the coulombic terms (first
transformation) followed by another transformation to relax the vdw terms
(in a separate transformation).
2) Then transform your CG1 and CG2 atoms to HB atoms in another set of
transformations.

During these transformation you have to choose the spacing between the
lambda values wisely based on some previous studies reported in literature
for similar type of transformation. Also the soft core parameters during
the (1) and (2) transformation should be set accrodingly.




On Thu, Jan 23, 2014 at 2:44 PM, Francesca Vitalini <
francesca.vitalini11 at gmail.com> wrote:

> Hi,
>
> thank you for the response,
>
>
> Il giorno 23/gen/2014, alle ore 06:35, bipin singh <bipinelmat at gmail.com>
> ha scritto:
>
> >> I am trying to perform a series of simulations where I gradually
> transform
> >> a valine capped amino acid (Ac-VAL-NHMe) into an alanine capped amino
> acid
> >> (Ac-ALA-NHMe), ie I have to turn on the HG and transform CG into HB. I
> have
> >> defined 6 lambda values at which I would like to equilibrate my system
> >> (lambda=0,0.2,0.40.6,0.8,1)
> >
> > Would it not be better if you turn off the vdw and coul. terms of  CG (1
> > atom) and HG (two atoms) and transform one of the HG to HB while doing
> the
> > transformation from valine to alanine, rather than mutating the CG to HB.
> >
>
> If I transform a HG into a HB then won't I have problems with the bonds
> definitions? How should I change my topology in order to define that in the
> end state the mutated HG has to be bound to CB instead?
>
> Thank you very much.
>
> Francesca
>
> >
> > On Wed, Jan 22, 2014 at 9:40 PM, Francesca Vitalini <
> > francesca.vitalini11 at gmail.com> wrote:
> >
> >> Dear Gromacs users,
> >>
> >> I am trying to perform a series of simulations where I gradually
> transform
> >> a valine capped amino acid (Ac-VAL-NHMe) into an alanine capped amino
> acid
> >> (Ac-ALA-NHMe), ie I have to turn on the HG and transform CG into HB. I
> have
> >> defined 6 lambda values at which I would like to equilibrate my system
> >> (lambda=0,0.2,0.40.6,0.8,1)
> >>
> >> I am using GROMACS.4.5.5. and OPLS-AA force field.
> >>
> >> I have been following the tutorial by Justin Lemkul
> >>
> >>
> http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin/gmx-tutorials/free_energy/index.html
> >>
> >> I have problems with the nvt equilibration with lambda=0.8 and lambda=1
> as
> >> my systems systematically crashes. I believe that the reason is in the
> >> tuning of the parameters sc-alpha sc-power sc-sigma in the mdp file.
> >> I understand from the manual that sc_sigma should not be the default
> option
> >> as I have 6 hydrogens disappearing, and consequently the sc-alpha
> sc-power
> >> parameters will be effected as the interactions of the hydrogens in the
> >> soft core state will get stronger.
> >> However I have no experience in how to select these parameters and would
> >> like to ask, if anyone has more experience, which could be plausible
> >> values.
> >>
> >> Thank you very much,
> >>
> >> Francesca
> >>
> >>
> >> For completeness I attach here a copy of my topology and mdp files.
> >>
> >>
> >> TOPOLOGY
> >>
> >>
> >> ; Include forcefield parameters
> >> #include "oplsaa.ff/forcefield.itp"
> >>
> >> [ moleculetype ]
> >> ; Name            nrexcl
> >> Protein_chain_A     3
> >>
> >> [ atoms ]
> >> ;   nr       type  resnr residue  atom   cgnr     charge       mass
> >> typeB    chargeB      massB
> >> ; residue   1 ACE rtp ACE  q  0.0
> >>     1   opls_135      1    ACE    CH3      1      -0.18     12.011
> >> opls_135     -0.18     12.011
> >>     2   opls_140      1    ACE   HH31      1       0.06      1.008
> >> opls_140      0.06      1.008
> >>     3   opls_140      1    ACE   HH32      1       0.06      1.008
> >> opls_140      0.06      1.008
> >>     4   opls_140      1    ACE   HH33      1       0.06      1.008
> >> opls_140      0.06      1.008
> >>     5   opls_235      1    ACE      C      2        0.5     12.011
> >> opls_235       0.5     12.011
> >>     6   opls_236      1    ACE      O      2       -0.5    15.9994
> >> opls_236      -0.5    15.9994
> >> ; residue   2 VAL rtp VAL  q  0.0
> >>     7   opls_238      2    VAL      N      3       -0.5    14.0067
> >> opls_238      -0.5    14.0067
> >>     8   opls_241      2    VAL      H      3        0.3      1.008
> >> opls_241       0.3      1.008
> >>     9  opls_224B      2    VAL     CA      3       0.14     12.011
> >> opls_224B      0.14     12.011
> >>    10   opls_140      2    VAL     HA      3       0.06      1.008
> >> opls_140      0.06      1.008
> >>    11   opls_137      2    VAL     CB      4      -0.06     12.011
> >> opls_137     -0.18     12.011
> >>    12   opls_140      2    VAL     HB      4       0.06      1.008
> >> opls_140      0.06      1.008
> >>    13   opls_135      2    VAL    CG1      5      -0.18     12.011
> >> opls_140      0.06      1.008
> >>    14   opls_140      2    VAL   HG11      5       0.06      1.008
> >> opls_140      0.00      0.000
> >>    15   opls_140      2    VAL   HG12      5       0.06      1.008
> >> opls_140      0.00      0.000
> >>    16   opls_140      2    VAL   HG13      5       0.06      1.008
> >> opls_140      0.00      0.000
> >>    17   opls_135      2    VAL    CG2      6      -0.18     12.011
> >> opls_140      0.06      1.008
> >>    18   opls_140      2    VAL   HG21      6       0.06      1.008
> >> opls_140      0.00      0.000
> >>    19   opls_140      2    VAL   HG22      6       0.06      1.008
> >> opls_140      0.00      0.000
> >>    20   opls_140      2    VAL   HG23      6       0.06      1.008
> >> opls_140      0.00      0.000
> >>    21   opls_235      2    VAL      C      7        0.5     12.011
> >> opls_235       0.5     12.011
> >>    22   opls_236      2    VAL      O      7       -0.5    15.9994
> >> opls_236      -0.5    15.9994
> >> ; residue   3 NAC rtp NAC  q  0.0
> >>    23   opls_238      3    NAC      N      8       -0.5    14.0067
> >> opls_238      -0.5    14.0067
> >>    24   opls_241      3    NAC      H      8        0.3      1.008
> >> opls_241       0.3      1.008
> >>    25   opls_242      3    NAC    CH3      9       0.02     12.011
> >> opls_242      0.02     12.011
> >>    26   opls_140      3    NAC   HH31      9       0.06      1.008
> >> opls_140      0.06      1.008
> >>    27   opls_140      3    NAC   HH32      9       0.06      1.008
> >> opls_140      0.06      1.008
> >>    28   opls_140      3    NAC   HH33      9       0.06      1.008
> >> opls_140      0.06      1.008
> >>
> >> [ bonds ]
> >> ;  ai    aj funct            c0            c1            c2
>  c3
> >>    1     2     1
> >>    1     3     1
> >>    1     4     1
> >>    1     5     1
> >>    5     6     1
> >>    5     7     1
> >>    7     8     1
> >>    7     9     1
> >>    9    10     1
> >>    9    11     1
> >>    9    21     1
> >>   11    12     1
> >>   11    13     1
> >>   11    17     1
> >>   13    14     1
> >>   13    15     1
> >>   13    16     1
> >>   17    18     1
> >>   17    19     1
> >>   17    20     1
> >>   21    22     1
> >>   21    23     1
> >>   23    24     1
> >>   23    25     1
> >>   25    26     1
> >>   25    27     1
> >>   25    28     1
> >>
> >> [ pairs ]
> >> ;  ai    aj funct            c0            c1            c2
>  c3
> >>    1     8     1
> >>    1     9     1
> >>    2     6     1
> >>    2     7     1
> >>    3     6     1
> >>    3     7     1
> >>    4     6     1
> >>    4     7     1
> >>    5    10     1
> >>    5    11     1
> >>    5    21     1
> >>    6     8     1
> >>    6     9     1
> >>    7    12     1
> >>    7    13     1
> >>    7    17     1
> >>    7    22     1
> >>    7    23     1
> >>    8    10     1
> >>    8    11     1
> >>    8    21     1
> >>    9    14     1
> >>    9    15     1
> >>    9    16     1
> >>    9    18     1
> >>    9    19     1
> >>    9    20     1
> >>    9    24     1
> >>    9    25     1
> >>   10    12     1
> >>   10    13     1
> >>   10    17     1
> >>   10    22     1
> >>   10    23     1
> >>   11    22     1
> >>   11    23     1
> >>   12    14     1
> >>   12    15     1
> >>   12    16     1
> >>   12    18     1
> >>   12    19     1
> >>   12    20     1
> >>   12    21     1
> >>   13    18     1
> >>   13    19     1
> >>   13    20     1
> >>   13    21     1
> >>   14    17     1
> >>   15    17     1
> >>   16    17     1
> >>   17    21     1
> >>   21    26     1
> >>   21    27     1
> >>   21    28     1
> >>   22    24     1
> >>   22    25     1
> >>   24    26     1
> >>   24    27     1
> >>   24    28     1
> >>
> >> [ angles ]
> >> ;  ai    aj    ak funct            c0            c1
> >> c2            c3
> >>    2     1     3     1
> >>    2     1     4     1
> >>    2     1     5     1
> >>    3     1     4     1
> >>    3     1     5     1
> >>    4     1     5     1
> >>    1     5     6     1
> >>    1     5     7     1
> >>    6     5     7     1
> >>    5     7     8     1
> >>    5     7     9     1
> >>    8     7     9     1
> >>    7     9    10     1
> >>    7     9    11     1
> >>    7     9    21     1
> >>   10     9    11     1
> >>   10     9    21     1
> >>   11     9    21     1
> >>    9    11    12     1
> >>    9    11    13     1
> >>    9    11    17     1
> >>   12    11    13     1
> >>   12    11    17     1
> >>   13    11    17     1
> >>   11    13    14     1
> >>   11    13    15     1
> >>   11    13    16     1
> >>   14    13    15     1
> >>   14    13    16     1
> >>   15    13    16     1
> >>   11    17    18     1
> >>   11    17    19     1
> >>   11    17    20     1
> >>   18    17    19     1
> >>   18    17    20     1
> >>   19    17    20     1
> >>    9    21    22     1
> >>    9    21    23     1
> >>   22    21    23     1
> >>   21    23    24     1
> >>   21    23    25     1
> >>   24    23    25     1
> >>   23    25    26     1
> >>   23    25    27     1
> >>   23    25    28     1
> >>   26    25    27     1
> >>   26    25    28     1
> >>   27    25    28     1
> >>
> >> [ dihedrals ]
> >> ;  ai    aj    ak    al funct            c0            c1
> >> c2            c3            c4            c5
> >>    2     1     5     6     3
> >>    2     1     5     7     3
> >>    3     1     5     6     3
> >>    3     1     5     7     3
> >>    4     1     5     6     3
> >>    4     1     5     7     3
> >>    1     5     7     8     3
> >>    1     5     7     9     3
> >>    6     5     7     8     3
> >>    6     5     7     9     3
> >>    5     7     9    10     3
> >>    5     7     9    11     3
> >>    5     7     9    21     3
> >>    8     7     9    10     3
> >>    8     7     9    11     3
> >>    8     7     9    21     3
> >>    7     9    11    13     3    dih_VAL_chi1_N_C_C_C
> >>    7     9    11    17     3    dih_VAL_chi1_N_C_C_C
> >>   21     9    11    13     3    dih_VAL_chi1_C_C_C_CO
> >>   21     9    11    17     3    dih_VAL_chi1_C_C_C_CO
> >>    7     9    11    12     3
> >>   10     9    11    12     3
> >>   10     9    11    13     3
> >>   10     9    11    17     3
> >>   21     9    11    12     3
> >>    7     9    21    22     3
> >>    7     9    21    23     3
> >>   10     9    21    22     3
> >>   10     9    21    23     3
> >>   11     9    21    22     3
> >>   11     9    21    23     3
> >>    9    11    13    14     3
> >>    9    11    13    15     3
> >>    9    11    13    16     3
> >>   12    11    13    14     3
> >>   12    11    13    15     3
> >>   12    11    13    16     3
> >>   17    11    13    14     3
> >>   17    11    13    15     3
> >>   17    11    13    16     3
> >>    9    11    17    18     3
> >>    9    11    17    19     3
> >>    9    11    17    20     3
> >>   12    11    17    18     3
> >>   12    11    17    19     3
> >>   12    11    17    20     3
> >>   13    11    17    18     3
> >>   13    11    17    19     3
> >>   13    11    17    20     3
> >>    9    21    23    24     3
> >>    9    21    23    25     3
> >>   22    21    23    24     3
> >>   22    21    23    25     3
> >>   21    23    25    26     3
> >>   21    23    25    27     3
> >>   21    23    25    28     3
> >>   24    23    25    26     3
> >>   24    23    25    27     3
> >>   24    23    25    28     3
> >>
> >> [ dihedrals ]
> >> ;  ai    aj    ak    al funct            c0            c1
> >> c2            c3
> >>    1     7     5     6     1    improper_O_C_X_Y
> >>    5     9     7     8     1    improper_Z_N_X_Y
> >>    9    23    21    22     1    improper_O_C_X_Y
> >>   21    25    23    24     1    improper_Z_N_X_Y
> >>
> >> ; Include Position restraint file
> >> #ifdef POSRES
> >> #include "Ac_V_NHMe/Ac_V_NHMe.itp"
> >> #endif
> >>
> >> ; Include water topology
> >> #include "oplsaa.ff/tip3p.itp"
> >>
> >> #ifdef POSRES_WATER
> >> ; Position restraint for each water oxygen
> >> [ position_restraints ]
> >> ;  i funct       fcx        fcy        fcz
> >>   1    1       1000       1000       1000
> >> #endif
> >>
> >> ; Include topology for ions
> >> #include "oplsaa.ff/ions.itp"
> >>
> >> [ system ]
> >> ; Name
> >> ACETYL-VALINE-METHYLAMIDE in water
> >>
> >> [ molecules ]
> >> ; Compound        #mols
> >> Protein_chain_A     1
> >> SOL               682
> >>
> >>
> >>
> >>
> >> NVT_VDW-Q.mdp
> >> here I intend to switch off the coulomb interactions.
> >>
> >> ; Run control
> >> integrator               = sd       ; Langevin dynamics
> >> tinit                    = 0
> >> dt                       = 0.002
> >> nsteps                   = 10000    ; 20 ps
> >> nstcomm                  = 100
> >> ; Output control
> >> nstxout                  = 10000
> >> nstvout                  = 10000
> >> nstfout                  = 0
> >> nstlog                   = 1000
> >> nstenergy                = 500
> >> nstxtcout                = 500
> >> xtc-precision            = 500
> >> xtc_grps                 = Protein
> >> ; Neighborsearching and short-range nonbonded interactions
> >> nstlist                  = 5
> >> ns_type                  = grid
> >> pbc                      = xyz
> >> rlist                    = 1.0
> >> ; Electrostatics
> >> coulombtype              = PME
> >> rcoulomb                 = 1.0
> >> ; van der Waals
> >> vdw-type                 = switch
> >> rvdw-switch              = 0.8
> >> rvdw                     = 0.9
> >> ; Apply long range dispersion corrections for Energy and Pressure
> >> DispCorr                  = EnerPres
> >> ; Spacing for the PME/PPPM FFT grid
> >> fourierspacing           = 0.12
> >> ; EWALD/PME/PPPM parameters
> >> pme_order                = 6
> >> ewald_rtol               = 1e-06
> >> epsilon_surface          = 0
> >> optimize_fft             = no
> >> ; Temperature coupling
> >> ; tcoupl is implicitly handled by the sd integrator
> >> tc_grps                  = system
> >> tau_t                    = 1.0
> >> ref_t                    = 300
> >> ; Pressure coupling is off for NVT
> >> Pcoupl                   = No
> >> tau_p                    = 0.5
> >> compressibility          = 4.5e-05
> >> ref_p                    = 1.0
> >> ; Free energy control stuff
> >> free_energy              = yes
> >> init_lambda              = 0.8
> >> delta_lambda             = 0
> >> sc-alpha                 = 0.5
> >> sc-power                 = 1.0
> >> sc-sigma                 = 0.005
> >> couple-moltype           = Protein_chain_A      ; name of moleculetype
> to
> >> decouple
> >> couple-lambda0           = vdw-q      ; only van der Waals interactions
> >> couple-lambda1           = vdw     ; turn off everything, in this case
> only
> >> vdW
> >> couple-intramol          = no
> >> ; Generate velocities to start
> >> gen_vel                  = yes
> >> gen_temp                 = 300
> >> gen_seed                 = -1
> >> ; options for bonds
> >> constraints              = h-bonds  ; we only have C-H bonds here
> >> ; Type of constraint algorithm
> >> constraint-algorithm     = lincs
> >> ; Do not constrain the starting configuration
> >> continuation             = no
> >> ; Highest order in the expansion of the constraint coupling matrix
> >> lincs-order              = 12
> >>
> >>
> >> NVT_VDW.mdp
> >> here I intend to switch off the Van der Waals interactions. The only
> >> difference with the previous file is in the Free energy parameters, so I
> >> will report only those.
> >>
> >> ; Free energy control stuff
> >> free_energy              = yes
> >> init_lambda              = 0.8
> >> delta_lambda             = 0
> >> sc-alpha                 = 0.5
> >> sc-power                 = 1.0
> >> sc-sigma                 = 0.3
> >> couple-moltype           = Protein_chain_A      ; name of moleculetype
> to
> >> decouple
> >> couple-lambda0           = vdw      ; only van der Waals interactions
> >> couple-lambda1           = none     ; turn off everything, in this case
> >> only vdW
> >> --
> >> Gromacs Users mailing list
> >>
> >> * Please search the archive at
> >> http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before
> >> posting!
> >>
> >> * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
> >>
> >> * For (un)subscribe requests visit
> >> https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or
> >> send a mail to gmx-users-request at gromacs.org.
> >>
> >
> >
> >
> > --
> >
> >
> >
> > *--------------------Thanks and Regards,Bipin Singh*
> > --
> > Gromacs Users mailing list
> >
> > * Please search the archive at
> http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before
> posting!
> >
> > * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
> >
> > * For (un)subscribe requests visit
> > https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or
> send a mail to gmx-users-request at gromacs.org.
>
> --
> Gromacs Users mailing list
>
> * Please search the archive at
> http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before
> posting!
>
> * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
>
> * For (un)subscribe requests visit
> https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or
> send a mail to gmx-users-request at gromacs.org.
>



-- 



*--------------------Thanks and Regards, Bipin Singh*


More information about the gromacs.org_gmx-users mailing list