From: Evgeny Bulat (jack.bulat_at_gmail.com)
Date: Tue Jun 05 2012 - 13:00:52 CDT
Hey everyone!
I'm pretty green on the VMD/NAMD front, and have an issue I can't seem to
resolve. For my work, I'm trying to construct a glycopeptide and model it
in NAMD. I used CHARMM and the CHARMM36 toppar package (containing the
_carb and _glycopeptide toppars, as well as the more established protein
CHARMM22 toppars) to construct it. Once that's done, I shift to VMD and use
the Psfgen 1.6 plugin to create segments using the glycopeptide PDB and
others, write psf/pdb files for the system, solvate, and ionize. Then I run
it out on NAMD. I turned to VMD for building my system because it's much
easier for me to work with than CHARMM.
However, I've run into several issues related specifically to the Psfgen
plugin. First of all, it models the glyco-adduct as just another residue.
In the CHARMM-originating PDB that I use to construct the segment in VMD,
the monosaccharide sugar is placed as the last (C-terminal) moiety. Psfgen
looks at this and attempts to connect the second-to-last moiety (a real
amino-acid residue, the "true" peptide C-terminus) to the sugar - even if,
structurally, the sugar is connected to a completely different residue and
not to that one. As a result, several things happen:
1) Upon segment construction, the automatic CTER patch is ignored since,
according to Psfgen, the sugar is at the C-terminus and that doesn't have
the CA atom.
2) There's a highly distorted amide bond that's made automatically placed
between the C-terminal amino-acid and the way-distant sugar adduct (which,
conveniently/inconveniently, carries an amino function).
As a way of tricking Psfgen to cooperate, I have resolved these issues by
doing the following (assuming that Psfgen is imported and all the
aforementioned par files are called):
segment GP { pdb glycopeptide_xtra.pdb } ;# devise a segment from a
glycopeptide carrying an extra amino-acid at the C-term (11 amino-acids
long + the sugar = 12 "residues")
patch SGPB GP:3 GP:12 ;# patch a glycosidic bond between amino-acid 3 and
the sugar (in the pdb, "amino-acid" 12)
delatom GP 11 ;#remove the extra amino-acid at the C-term, since that's
automatically connected to the sugar as well by Psfgen
patch CTER GP:10 ;#generate a new C-terminus on the new last amino-acid
coordpdb glycopeptide_normal.pdb ;#get coords from a PDB of a glycopeptide
of normal length (10 amino-acids long + the sugar = 11 "residues")
guesscoord
writepsf blah.psf
writepdb blah.pdb
I then open the molecule:
mol new blah.psf
mol addfile blah.pdb
And it looks fine! To reiterate, I had to truncate the C-terminal
amino-acid since Psfgen would otherwise put in an absurd amide bond between
the C-terminus and the sugar (which I definitely don't want). If you follow
the numbering, it would go from 1-10 and then jump to 12 for the new
structure (12 being the sugar). For glycopeptide_normal.pdb, I actually had
to go in and renumber the sugar as residue 12 rather than 11 (the default,
by simple order) for the above trick to work.
Things seemed to go well. However, a serious issue arose in NAMD. And not
for the second issue with the Psfgen plugin: the manual CTER patch doesn't
work right! At the bottom of this post is my NAMD configuration file for
running a system consisting of a solvated glycopeptide (for simplicity's
sake). The errors I repeatedly get are one or more of the following:
1) ERROR: Atom 4602 velocity is 100492 -42002.8 114437 (limit is 6000)
ERROR: Atoms moving too fast; simulation has become unstable.
ERROR: Constraint failure in RATTLE algorithm for atom 4849!
ERROR: Constraint failure; simulation has become unstable.
2) Signal: segmentation violation
3) FATAL ERROR: Periodic cell has become too small for original patch grid!
Possible solutions are to restart from a recent checkpoint,
increase margin, or disable useFlexibleCell for liquid simulation.
These issues consistently come up after minimization finishes and
equilibration kicks in. Error 3 is actually resolvable and not totally
applicable to my real system. However, errors 1 and 2 are. I suspect that
Psfgen's CTER patch is to blame because:
1) The atoms that are messed up in error 1 are ALWAYS either the CA or the
OT1/OT2 atoms of the C-terminal amino-acid (without exception).
2) I made a segment of a regular peptide (without the sugar) with and
without manually applying the CTER patch after segment generation. When I
applied the patch manually, the same problems arose (so it's not the
glycosylation per se that's to blame). When I didn't apply the patch,
everything ran just fine! Psfgen autopatched the C-terminus during segment
generation, and that worked fine. However, it won't autopatch when there's
a sugar because that's treated at the C-terminal residue.
3) I managed to get a glimpse into the geometry of the manually patched
C-terminus BRIEFLY into equilibration of the peptide system. Those OT1/OT2
atoms for doing some crazy acrobatics! The proper trigonal planar geometry
was definitely not preserved by Psfgen!
I have tried the following already: 1) changing the order of the sugar in
the PDB; 2) imposing a step-wise temperature increase for my NAMD
simulation; 3) starting with equilibration right away (no minimizing).
The little bit of success that I saw is in running my CHARMM-created
PSF/PDB combo alone in NAMD (after converting the PSF from CHARMM to
X-Plor, of course). It ran fine, the geometry was correct, beautiful.
However:
1) Running a glycopeptide solely by itself is obviously not exactly what I
want to do.
2) VMD doesn't recognize the CHARMM-originating PSF for the purposes of
solvation and writing new PSFs based on global atom selection. Alongside,
Psfgen can't read any of my CHARMM peptide and glycopeptide PSF when I try
to use readpsf.
At this point, I'm a little lost. I hate to admit it, but I might have to
create the ENTIRE system within CHARMM! To avoid this, I welcome any advice
and suggestions as to what to try next. BTW, sorry for the length, but I
just wanted to spell everything out. Thank you very much in advance!
##############################
###############################
## JOB DESCRIPTION ##
#############################################################
# Equilibration after minimization was performed
#############################################################
## ADJUSTABLE PARAMETERS ##
#############################################################
structure gp_solvated.psf
coordinates gp_solvated.pdb
set temperature 310
set outputname gp_solvated_eq
firsttimestep 0
# restart conditions
if {0} {
set inputname -
binCoordinates -
binVelocities -
extendedSystem -
}
#############################################################
## SIMULATION PARAMETERS ##
#############################################################
# Input
paraTypeCharmm on
parameters par_all27_prot_lipid_na.inp
parameters par_all36_carb_namd.prm
parameters toppar_all36_carb_model.str
parameters toppar_all36_glycopeptide.str
temperature $temperature
# Force-Field Parameters
exclude scaled1-4
1-4scaling 1.0
cutoff 12.0
switching on
switchdist 10.0
pairlistdist 13.5
# Integrator Parameters
timestep 2.0 ;# 2fs/step
rigidBonds all ;# needed for 2fs steps
nonbondedFreq 1
fullElectFrequency 2
stepspercycle 10
# Constant Temperature Control
langevin on ;# do langevin dynamics
langevinDamping 5 ;# damping coefficient (gamma) of 5/ps
langevinTemp $temperature
langevinHydrogen off ;# don't couple langevin bath to hydrogens
if {1} {
# Periodic Boundary Conditions
cellBasisVector1 65.0 0. 0.0
cellBasisVector2 0.0 51.0 0.0
cellBasisVector3 0.0 0 39.0
cellOrigin 0.000006836525699327467 -0.000011816733604064211
-0.00000097775409813039
}
wrapAll on
# PME (for full-system periodic electrostatics)
PME yes
PMEGridSpacing 1.0
#manual grid definition
#PMEGridSizeX 150
#PMEGridSizeY 150
#PMEGridSizeZ 128
# Constant Pressure Control (variable volume)
if {1} {
useGroupPressure yes ;# needed for rigidBonds; specifies whether ints
involving Hs should be counted for all Hs or b/w groups
useFlexibleCell no
useConstantArea no
langevinPiston on
langevinPistonTarget 1.01325 ;# in bar -> 1 atm
langevinPistonPeriod 100.0
langevinPistonDecay 50.0
langevinPistonTemp $temperature
}
# Output
outputName $outputname
restartfreq 5000 ;# 1000steps = every 1ps
dcdfreq 100
xstFreq 100
outputEnergies 100
outputPressure 100
#############################################################
## EXTRA PARAMETERS ##
#############################################################
#Fixed atom constraints
if {0} {
fixedAtoms on
fixedAtomsFile
fixedAtomsCol B
fixedAtomsForces on
}
#############################################################
## EXECUTION SCRIPT ##
#############################################################
# Minimization
minimize 10000
reinitvels $temperature
run 25000 ;# 50ps
#############################################################
## END #
#############################################################
-Jack Bulat
Depts. of Molecular and Cell Biology and Bioengineering, UC Berkeley
This archive was generated by hypermail 2.1.6 : Tue Dec 31 2013 - 23:22:04 CST