From:	Jeffrey.Nauss()at()UC.Edu (Jeffrey L. Nauss)
Date:	Fri, 17 Jan 1997 10:58:13 -0500
Subject:	Summary: Removing charges for an MD simulation

Recently, I sent out the following post:

+++++++++++ BEGIN ORIGINAL MESSAGE ++++++++++++++++

One of the problems one has in MD simulations is that long runs must be
performed to see anything significant, for example, and in my specific case,
determination of NMR order parameters or conformational searches.  Performing
the run in explicit solvent with periodic boundary conditions is undoubtably
the best way to perform the simulations.  However, these simulations require
large amounts of CPU time and disk space for the large history files.

To get around these problems, people have run simulations in vacuo.  This was
done in the early days of MD simulations but with parameter sets that often
took into account the absence of solvent.  Furthermore, electrostatic screening
was accomplished by using a distance dependent dielectric constant.  Still, a
continuing problem that I have found during an in vacuo run is the collapse of
the molecule.  Proteins will shrink (in terms of the radius of gyration) and
side chains will fall onto the surface of the protein.  Peptides often fold up
into themselves.  To some extent, this collapse is a result of coulmbic
interactions that are not properly shielded.

Ornstein's group tried to get around this problem by rendering acidic, basic,
N-termini, and C-termini functional groups net neutral (JBSD vol 9, page 935
(1992)).  It appeared to be somewhat successful.

Has anyone tried MD simulations, particularly with small molecules like
peptides and carbohydrates, with a similar approach?   I am thinking perhaps of
performing the simulation with the Coulombic interactions scaled down or even
shut off.  Of course, one must assume that the Coulombic interactions will not
play a major role in the dynamics.  How valid is that assumption will depend on
the system studied, I would think.

Comments will be appreciated.  References will be even better.

+++++++++++ END ORIGINAL MESSAGE ++++++++++++++++

There seems to be a fair of amount of thought on the subject.  And there
seems to be some good results with various techniques.  At least enough
for me to give it a try.

Here is a summary of the responses.  Thank you to all who responded.

========================================================================

From: bear[ AT ]ellington.pharm.arizona.edu (Soaring Bear)

	Are you aware of the approach of using higher
dielectric r?   There's been a bit of literature on
using 4 instead of 1.   I could dig up some refs if
you want.

>Forgot to mention that.  Do you have a ref or two for it?

comparison in:
J Biomolec Struct & Dynam 11:429-41 1993

nice discussion by Daggett/Kollman:
Biopolymers 31:285-304 1991
Chemica Scripta 29a:205-15 1989

mathematic treatment by Ramachandran:
Indian J Biochem 7:95-7 1970

========================================================================

From: D.van.der.Spoel (+ at +) chem.rug.NL

A conservative answer, just use the solvent and save yourself a lot of
trouble justifying simulations without solvent. You can do nanosecond
simulations of a protein in water on a reasonable workstation, or
even a Pentium Pro PC. If performance is crucial, check out the
GROMACS software, it is probably the fastest MD code on workstations.
http://rugmd0.chem.rug.nl/~gmx

========================================================================

From: Peter Grootenhuis 

I fully appreciate the problem you described with long MD runs.
It inspired us to generate a CHARMm-based force field for carbohydrates
(CHEAT95) in which the effect of solvent is implicitly included.
We also performed several MD-runs and typically got good results.
The references are: Grootenhuis & Haasnoot, Mol.Simul. 10 (1993) 75-95
and Kouwijzer & Grootenhuis, J. Phys. Chem 99 (1995) 13426-13436
The CHEAT95 force field is available at:
http://www.msi.com/support/quanta/cheat95.html

========================================================================

From: John Liebeschuetz 

	We've, in the past, carried out a number of in vacuo simulations of
drug/DNA complexes. No counterions were included so we scaled the phosphate
charges to -0.25 to simulate counterion shielding. This figure was rather ad
hoc but seemed to give reasonable results. Some of the ligands examined were
cationic. These were scaled to a charge of +0.5 over the whole molecule and
again seemed to give good results. In most of these simulations we still had to
restrain the DNA to stop it collapsing.
	The charges on the cationic ligand were originally ESP calculated in
AM1. We could have just halved the charges on each atom but we wanted to
maintain the charge differentials between electropositive and electronegative
atoms. This was done by carrying out the AM1 calculation on the unprotonated
and protonated species and then averaging the charges on each atom. The charge
on the extra proton was just halved.

> Have you a reference for this work?

Not for cases where the charges were partially removed. I do have one however
for a case where the DNA phosphates were made neutral.

N.L.Fregeau, Y. Wang, R.T.Pon, W.A.Wylie & J.W.Lown, J. Amer. Chem. Soc.,
117 (1995), 8917-8925.

> Did you manipulate the partial charges on the atoms to arrive at this scaled
> down value or did you scale all Coulombic interactions?

	for DNA, the charges on the DNA phosphate oxygens only were scaled
down. The charge on each atom was reduced by a quantity proportional to the
original negative charge. This method was used whether the group was to be
partially or fully neutralised.

 	The ligand was treated differently. When the ligand cationic charge was
scaled according to the method I described, the net effect was to alter atomic
charges in the vicinity of the cation centre only.

========================================================================

From: Pieter Stouten 

>Ornstein's group tried to get around this problem by rendering acidic,
>basic, N-termini, and C-termini functional groups net neutral (JBSD
>vol 9, page 935 (1992)). It appeared to be somewhat successful.
>
GROMOS has had a force field with net neutral groups for vacuum simulations
at least since 1984 (when I started using it), but probably much longer. It
definitely works better than full-blown charges. See e.g.:

  D.M.F. van Aalten, A. Amadei, R.P. Bywater, J.B.C. Findlay, H.J.C.
  Berendsen, C. Sander & P.F.W. Stouten, "A Comparison of Structural and
  Dynamic Properties of Different Simulation Methods Applied to SH3,"
  Biophys. J. 70 (1996) 684-692.

We also found that stochastic dynamics (SD, with friction and random
"kicks") worked better than straight MD. We also tried a simplistic
continuum solvation term described in:

  P.F.W. Stouten, C. Froemmel, H. Nakamura & C. Sander, "An Effective
  Solvation Term Based on Atomic Occupancies for Use in Protein
  Simulations," Mol. Simulation 10 (1993) 97-120.

That method is computationally cheap (20-50% increase in CPU time with
respect to vacuum simulations), did work well in BPTI simulations (in fact,
when I redid the simulations and was more careful in generating the
starting configuration the results were better even than published in the
paper), but not so well in the SH3 simulations.

========================================================================

From: Mark D Shenderovich 

I have discussed similar problems with fellow molecular modelers working
with peptides and proteins. I don't know a reference that compares
different electrostatic approximations in a systematic way, but I can
share with you results of our discussions and my own experience.

Charged groups of peptides in water are hydrated and interact through a
medium of a high dielectrics.  Using charged states of terminal and side
chain groups in vacuo tremendously overestimates their electrostatic
interactions and causes the collapse you've described. The charged groups
tend to form salt bridges or multi-H-bond networks. To my experience, a
distance dependent dielectrics makes the situation even worse since it
stronger punishes separation of opposite charges. You may try a
continuum solvation model which would increase CPU time by a reasonable
factor of 2-5 and will not require additional disk space. I don't know
which software you use, but several commercial packages include a
continuum solvation. Although most of these models are poorly parametrized
for charged species, it is better than nothing. Second, it is always safer
to use neutral form of ionizable groups, both in vacuo and with an
"envelope" hydration models (i.e. those which do not calculate effective
dielectric constants).   If you prefer simulations in vacuo, use a higher
dielectric constant (4 to 10) with neutral groups, or a distance dependent
dielectric with a coefficient around 4.0. I cannot prove these numbers,
they came from my and my colleagues' experience. Also, if possible with
your software, use aliphatic hydrogen parameters for neutral COOH groups.
This is recommended by people from Scheraga's group to avoid excessive
H-bonding of carboxyl. Finally, for conformational searches related to NMR
data I often use a continuum dielectric constant of the solvent (80 for
water, 45 for DMSO) with charged ionizable groups. Of course, this
reduces electrostatics for neutral groups inside the peptide/protein, but
they ether don't play a major role in intramolecular interactions or not
vary much from one conformation to another.

========================================================================

From: JAQ[ AT ]XRAY.BMC.UU.SE

I saw you ccl-posting. You might find the references J. Mol. Biol. (1985)
183: 461 and Biopolymers (1990) 30: 205 useful. We looked at some of
the issues you mention there, and there may be some other useful references
for you as well.

========================================================================

From: "Thomas M O'Connell" 

We recently published a study on the helix-coil transition
for short peptides that may be of interest to you.  In this
study we ran dynamics/free energy simulations in vacuo with
the dielectric constant set to 5, 25 and infinity (i.e. all
charges turned off).  This of course had definite effects on
the thermodynamics of the peptide conformation and folding.
The reference is : Wang, O'Connell, Tropsha and Hermans
                   Biopolymers  39, 479 (1996)

Also of interest is a study that came out a while ago from
Bernie Brooks in which he ran some dynamics on BPTI (if I remember
correctly) in which all charges were removed.  I can't seem to
find this in my files so you'll have to do a search, but not
surprizingly the folded system was stable.

========================================================================

From: case - at - riscsm.scripps.EDU (David Case)

We've thought a lot about doing this, but for some reason have never actually
tried it....For peptides, some models of solvation effects almost perfectly
cancel out the conformational dependence of the electrostatic term, i.e.
the sum of the (gas-phase) electrostatic term and the solvation term is
nearly constant as a function of conformations sampled during a solvated MD
run.  So it's certainly not a crazy idea.

No references yet, but we have a couple of papers submitted on this subject.

========================================================================

--
						Jeff Nauss

***********************************************************************
*  UU    UU             Jeffrey L. Nauss, PhD                         *
*  UU    UU             Director, Molecular Modeling Services         *
*  UU    UU             Department of Chemistry                       *
*  UU    UU CCCCCCC     University of Cincinnati                      *
*   UU  UU CCCCCCCC     Cincinnati, OH 45221-0172                     *
*    UUUU CC                                                          *
*         CC            Telephone: 513-556-0148    Fax: 513-556-9239  *
*         CC                                                          *
*          CCCCCCCC     e-mail: Jeffrey.Nauss -AatT- UC.Edu                  *
*           CCCCCCC     URL  http://www.che.uc.edu/~nauss
*
***********************************************************************



Similar Messages
01/09/1997:  Removing charges for an MD simulation
07/05/1997:  MD in various solvents. The Fortunate summary
08/01/1996:  Re: CCL:M:Heat of formation calculation using MOPAC.
09/01/1995:  Re:  CCL:M:Questions re dielec. const. & conformers
07/11/1994:  SUMMARY; MD & MC; Modelling of counterions effects ...
07/07/1995:  Summary: Opinions on the quantum charges
10/13/1992:  On the use of Cut-off Schemes in MD
04/04/1995:  Ewald sums
06/10/1993:  Solvation Models -- Time and Quality
11/26/1997:  tRNA modelling: summary


Raw Message Text