AMBER Archive (2009)Subject: Re: [AMBER] glycopeptide Parametrization
From: Karl Kirschner (kkirsch_at_scai.fraunhofer.de)
Date: Fri Jun 12 2009 - 10:06:17 CDT
Hi Carlo,
At a quick glance, the leap sequence of commands you have shown below (in
your most recent email) for linking the carbohydrate to the protein look
right to me.
In terms of the partial atomic charges, I would use the default charges that
come with the protein and carbohydrate residues, which are by default used
when you source their respective leaprc files. They have been determined in a
slightly different manner, which can be understood in their respective
publications.
We recommend that the default scaling factors (SCNB=2.0, SCEE=1.2) are used
when doing a simulation involving a protein and carbohydrates. However, this
has not been extensively tested and you should analyze your results to insure
that nothing "extrodinary" happened. (Based on what we have seen so far, the
use of the default scaling factors impacts the carbohydrate paraemeters by
lowering the rotational minima. Some small shifts in the relative minima are
also seen for the rotation about bonds). Right now, as far as I know, only
one set of scaling factors is used for your entire molecular system.
As for "overstabilization of the polar group-charged group interactions in
GLYCAM," I am not aware of any overstabilization. But that is something that
needs to be evaluated based on the simulation and comparison to experiment.
There is always a possiblity for this to occur with any force field and their
mixing, and Glycam is not immune to these issues. It could be that the
charges on the carbohydrates versus the protein are not 100% compatible due
to their slightly different development. If a formal charge is present, it
might manifest itself more. I would say that this is an area that could be
explored for a better understanding on how compatible they are with each
other. Hope that helps.
Cheers,
Karl
On Friday 12 June 2009 16:06, guardiani_at_fi.infn.it wrote:
> Dear Amber users,
>
> I would like to thank you all for your help and in particular
> I would like to thank Francois Dupradeau. Your idea to simplify
> the parametrization using only the Amber (ff99SB) force field
> and the GLYCAM force field is quite convincing.
>
> Now my parametrization is as simple as that (I am using the bondi
> radii because I want to try an implicit-solvent simulation with
> igb=7, but of course I will also run the explicit-solvent one):
>
> xleap -s -f $AMBERHOME/dat/leap/cmd/leaprc.ff99SB
>
> set default PBradii bondi
> source leaprc.glycam04
> x = sequence { NVAL GLU ARG NLN GLY HIS CSER }
> set x tail x.4.11
> set 1GB head 1GB.1.1
> y = sequence { x 1GB }
> savepdb y y.pdb
> saveamberparm y y.top y.crd
> quit
>
>
> Please let me know if this is what you meant. Which charges I will be
> using with this parametrization ? Charges based on the Connolly surface
> for the amino-acid residues and charges based on the CHELPG algorithm
> for the sugar ?
>
> As for the 1-4 scaling factors, I performed my previous simulations
> using the default values, e.g.
>
> SCNB = 2.0
> SCEF = 1.2
>
> Are these values appropriate for a glycopeptide simulation ?
> If I explicitly set these values in the .in input files of
> my REMD simulation, will the apply to both the peptide and
> the sugar parts ?
>
> Thank you very much again for your help.
>
> Best regards,
>
> Carlo
>
> Quoting FyD <fyd_at_q4md-forcefieldtools.org>:
> > 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
|