AMBER Archive (2009)

Subject: Re: [AMBER] glycopeptide Parametrization

From: FyD (fyd_at_q4md-forcefieldtools.org)
Date: Fri Jun 12 2009 - 02:05:51 CDT

Hi Carlo,

You provide a lot of information...

Here you mix a peptide part + a sugar part + am1-cc charges + GB +
.... These means you have an Amber force field for your protein (which
one ?) + a force field for the sugar (glycam) + am1-bcc charges (with
gaff ?)...

Moreover, which 1-4 scaling factors did you use for the non-bonding
interactions ? Amber FF and GLYCAM FF use different ones. Finally,
charges in Amber FF are based on the Connolly surface while GLYCAM FF
used the CHELPG algorithm.

I am not sure what to answer here. In these conditions, I would try to
"simplify" your computational conditions; i.e. using a minimum number
of different force fields; mixing ONLY Amber + Glycam FFs using
correct 1-4 scaling factors for the 1-4 interactions; This means
avoiding am1-bcc charges and using RESP charges based on conformations
you fully control. I would also use in a 1st step explicite solvent.

regards, Francois

> I am currently trying to perform simulations of a 21-residues long
> glycopeptide. The NMR structure of this peptide includes a beta-hairpin
> (residues 1-14) followed by a short alfa-helix. The orientation of the
> alfa-helix is rather variable, but in general it points away from the
> beta turn. On top of the beta turn there is an Asparagine linked to
> a beta-glucose unit through an N-glycosidic bond.
>
> Before describing my calculation procedure, I anticipate a bizarre
> result. The hairpin region produced by the REMD simulation in
> implicit solvent (igb=7), superposes reasonably well to the corrisponding
> NMR region. Conversely the region that was expected to adopt an helical
> structure, runs parallel to the C-terminal beta-strand. I noticed that
> in almost all conformations, this structure is stabilized by several
> H-bonds between the oxydrils of glucose and the C-terminal COO- group
> of the peptide. I therefore suspect that this kind of interactions may
> be overstabilized by Glycam. Do you have any experience of this ?
>
> First of all, here I show how I attained the topology of the glycopeptide.
>
> STEP 1
> -----------
>
> I built with xleap a new unit called ASG comprizing a capped
> asparagine linked
> to a beta-D_glucose molecule.
>
> xleap -s -f $AMBERHOME/dat/leap/cmd/leaprc.ff99
>
> source leaprc.glycam04
> set ACE tail ACE.1.5
> set ASN head ASN.1.1
> set ASN tail ASN.1.13
> set NME head NME.1.1
> x = sequence { ACE ASN NME }
> remove x x.2.11
> set x tail x.2.10
> set 1GB head 1GB.1.1
> y = sequence { x 1GB }
> savepdb y Asn_Glc.pdb
> quit
>
>
> STEP 2
> ------------
>
> As in the file Asn_Glc.pdb the two hydrogens bound to atom C6 of glucose are
> almost superposed, I manually corrected they coordinates. First of all
> I measured the correct bond distances C6-H16 and C6-H26 and the bond angle
> H16-C6-H26 in beta-D-glucose.
>
> xleap -s -f $AMBERHOME/dat/leap/cmd/leaprc.ff99
>
> source leaprc.glycam04
> x = copy 1GB
> savepdb x 1GB.pdb
> quit
>
> Using VMD I measured the parameters of interest:
>
> vmd 1GB.pdb
>
> measure bond {5 6} molid 0
> measure bond {5 7} molid 0
> measure angle {5 6 7} molid 0
>
> If we now define the coordinates of the two hydrogens as (x1,y1,z1) and
> (x2,y2,z2) and if we call theta the bond angle H16-C6-H26, and b the
> C-H bond length, it is possible to establish the correct relative position
> of the two hydrogens:
>
> x2 = x1 + 2*b*sin(theta/2)
> y2 = y1
> z2 = z1
>
> The edited ASG unit is saved in the file Asn_Glc_edit.pdb.
>
> STEP 3
> -----------
>
> Using VMD, the Asn_Glc_edit.pdb file is saved in mol2 format:
> Asn_Glc_edit.mol2
>
>
> STEP 4
> ------------
>
> Preparation of input file foe Gaussian run:
>
> antechamber -fi mol2 -fo gzmat -i Asn_Glc_edit.mol2 -o Asn_Glc_edit.pdb
>
> STEP 5
> ------------
>
> Gaussian run to optimize the structure and to compute the ESP charges.
> The input file for Gaussian has the form:
>
> --Link1--
> %Nproc=2
> %chk=molecule
> #HF/6-31G* SCF=Tight Test Pop=MK
> iop(6/33=2) iop(6/42=6) opt
> optimization of Asn_Glc
>
> HERE I PUT THE Z-MATRIX OF ASN_GLC
>
> The Gaussian job is started with the command:
>
> bsub -W 24:00 -n 2 -i Asn_Glc_edit.com -o Asn_Glc_edit.log
> -e Asn_Glc_edit.err g03
>
>
> STEP 6
> --------------
>
> Here we generate an .ac file where the ESP charges are converted to RESP
> charges, but the atomtypes are still those of the GAFF force field.
>
> antechamber -fi gout -fo ac -i Asn_Glc_edit.log
> -o Asn_Glc_edit.ac -c resp
>
>
> STEP 7
> ---------------
>
> Creation of a new .ac file where the GAFF atomtypes are converted to
> the AMBER atomtype:
>
> atomtype -i Asn_Glc_edit.ac -o Asn_Glc_edit_amber.ac -p amber
>
>
> STEP 8
> ---------------
>
> Preparation of a file mainchain.AsnGlc where in which we identify the atoms
> of the capping groups as atoms to be removed.
>
> HEAD_NAME N1
> TAIL_NAME C6
> MAIN_CHAIN C3
> OMIT_NAME C7
> OMIT_NAME H10
> OMIT_NAME H11
> OMIT_NAME H12
> OMIT_NAME N3
> OMIT_NAME H9
> OMIT_NAME C2
> OMIT_NAME O1
> OMIT_NAME C1
> OMIT_NAME H1
> OMIT_NAME H2
> OMIT_NAME H3
> PRE_HEAD_TYPE C
> POST_TAIL_TYPE N
> CHARGE 0.0
>
>
> STEP 9
> ---------------
>
> Preparation of prepi file.
>
> prepgen -i Asn_Glc_edit_amber.ac -o Asn_Glc_edit_amber.prepi
> -f prepi -m mainchain.AsnGlc -rn ASG -rf ASG.res
>
> STEP 10
> ---------------
>
> Preparation of frcmod file
>
> parmchk -i Asn_Glc_edit_amber.prepi -o Asn_Glc_edit_amber.frcmod -f prepi
>
>
> STEP 11
> ---------------
>
> Creation of topology file for implicit-solvent (igb=7) REMD simulation.
>
> xleap -s -f $AMBERHOME/dat/leap/cmd/leaprc.ff99SB
>
>
> set default PBradii bondi
> loadamberprep Asn_Glc_edit_amber.prepi
> loadamberparams Asn_Glc_edit_amber.frcmod
> x = sequence { NTHR PRO ARG VAL GLU ARG ASG GLY HIS SER VAL PHE LEU ALA
> PRO TYR GLY TRP MET VAL CLYS}
> saveamberparm x csf114_glc_ext.top csf114_glc_ext.crd
> quit
>
>
> =========================================================================
> REMD SIMULATION
> =========================================================================
>
> STEP 1: MINIMIZATION
> ----------------------
>
> Minimization
> &cntrl
> imin = 1,
> maxcyc = 400,
> ncyc = 200,
> ntb = 0,
> igb = 7,
> gbsa = 1,
> saltcon = 0.2,
> cut = 999.0,
> rgbmax=999.0,
> ntpr = 10,
> ntx = 1,
> /
>
> mpirun $sander -O -i minim.in -o minim.out -c csf114_glc_ext.crd -p
> csf114_glc_ext.top -r minim.rst > out
>
>
> STEP 2: EQUILIBRATION
> ------------------------
>
> We consider 20 different temperatures chosen in geometric progression
> T = 310, 324, 338, 353, 368, 384, 401, 418, 437, 456, 476, 497, 519,
> 541, 565, 590, 616, 642, 671, 700 K. Each replica is gradually heated
> to the corresponding temperature.
>
> Equilibration
> &cntrl
> imin = 0,
> cut = 999.0,
> rgbmax = 999.0,
> ig= 909,
> igb = 7,
> gbsa = 1,
> saltcon = 0.2,
> nstlim = 5000000,
> dt = 0.001,
> ntt = 3, gamma_ln = 1.0,
> tempi = 0.0, temp0 = 310.0,
> ntx = 1, irest = 0, ntb = 0,
> nscm = 50,
> ntpr = 100, ntwr = 100, ntwx = 100
> /
>
>
> STEP 3: REMD RUN
> --------------------
>
> We performed REMD 250 ns of REMD run proceding through restarts of 10 ns.
>
> MD run
> &cntrl
> imin = 0,
> cut = 999.0,
> rgbmax = 999.0,
> ig = 11377,
> igb = 7,
> gbsa = 1,
> saltcon = 0.2,
> nstlim = 2500,
> numexchg = 2000,
> dt = 0.001,
> ntt = 3, gamma_ln = 1.0,
> tempi = 310.0, temp0 = 310.0,
> ntx = 5, irest = 1, ntb = 0,
> nscm = 50,
> repcrd = 0,
> ntpr = 500, ntwr = 500, ntwx = 500
> /
>
>
> What did I do wrong ? Is there anything that I may try to
> improve the parametrization of the glycopeptide ? Is there
> any known overstabilization of the polar group-charged group
> interactions in GLYCAM ?
>
> Thank you very much for your help.
>
> Best regards,
>
> Carlo Guardiani
>

_______________________________________________
AMBER mailing list
AMBER_at_ambermd.org
http://lists.ambermd.org/mailman/listinfo/amber