Date: Thu, 21 May 1998 15:42:37 -0400
From: "David A. Pearlman"
Organization: Vertex Pharmaceuticals, Inc.
X-Mailer: Mozilla 3.02 (X11; I; IRIX 6.2 IP22)
MIME-Version: 1.0
To: amber@cgl.ucsf.EDU
Subject: Re: GIBBS alternative topology definitions?
Tom Pochapsky wrote:
>
> I am interested in using an alternative to the standard method of
> implementing parameter changes in GIBBS calculations as discussed
> in the GIBBS manual, section 1.8.7, to wit.....
>
> >1.8.7. VII. Changing parameters versus dual topologies
> >
> > In "standard" operation, free energy changes in GIBBS are
> > effected by transforming the potential representative of state
> > 1 to that representative of state 2. The topology of the system
> > does not change. To make atoms non-interacting at one of the
> > endpoints, they are assigned zeroed non-bond and electrostatic
> > parameters at this endpoint.
>
> > The improved mixing rules which can be used in GIBBS Ver-
> > sion 4 (IOLEPS = 0, line 10) allow a second method to be used.
> > One result of these new mixing rules is that if any pair of
> > atoms "exist" only at mutually exclusive endpoints (e.g. atom i
> > exists in state 1 but not state 2; atom j exists in state 2,
> > but not in state 1), then effectively no non-bonded interac-
> > tions are ever calculated between them. This means that, in
> > lieu of the "standard" method which uses a single topology, we
> > can define dual topologies, one corresponding to the lambda = 0
> > endpoint, and the other corresponding to the lambda = 1 end-
> > point. For the former topology, all non-bonded parameters would
> > be defined to be 0 in the lambda = 1 state. Similarly, all non-
> > bonded parameters for the latter topology would be 0 at lambda
> > = 0. The two topologies would then never "see" each other at
> > intermediate values of lambda. Defining dual topologies can
> > aid in performing free energy calculations where bond connec-
> > tivities must change. Dual topologies is the method incorpo-
> > rated into the "CHARMM" program.
>
> Has this actually been implemented? If so, how does one go about
> defining the second (lambda=0) topology? Any help on this would be
> appreciated.
>
Yes, you can implement a "Charmm" type change, if that's what you
want. To do so, you need to define topologies for both the beginning
and ending states in the prep file (or during Leap, if you're using
that to generate the topology file). As an example, if you wanted to
change Glycine to Alanine, you would define the topology change as:
COO- COO-
| |
H3N - C - H ---> H3N - C - H
/ \ / \
H DC(DH)3 DH CH3
Then, you would define the parameters such that:
1) DC and DH are "dummy" non-interacting atoms (epsilon = 0, and
charge = 0);
2) All internal coordinates that arise from having the beginning
and ending sidechains on the C-alpha have 0 force constants.
In this case, this would mean setting K=0 for the parameters
for the H-C-DC and DH-C-C valence angles made by the atoms of the
sidechain.
Now here's the tricky part: We can't just set the force constant
for the H-C-DC or DH-C-C valence angle to 0. Why? Because if
we define the atoms with standard atom types, not only is the
angle between the two sidechains defined by these types, the
angle between the sidechain and _other_ substituents of
C-alpha are also defined by these. And those angles should
retain standard (non 0) valence angle parameters. What this
means is that you need to define the atom types of the
atoms in the sidechain uniquely. In turn, this will require going
back to the parameter file and including all the parameters
that this will necessitate.
For example, if we define the sidechain H to be HX (instead of
the usual HC), then we need to define all the following parameters
which would not usually be in the parameters file:
HX-CT
HX-CT-DC
HX-CT-N
HX-CT-HC
HX-CT-C
HX-CT-N-H
HX-CT-C-O
HX-CT-C-OH
etc.
as well as
HX-CT-DC with K=0, as mentioned.
If you use the Amber (single topology) approach, the change is
defined as
COO- COO-
| |
H3N - C - H ---> H3N - C - H
| |
H(DH)3 CH3
Here you do not need to define a special hydrogen type because
there are no internal coordinates (e.g. valence angles) that need
to be zeroed out.
The bottom line is: The Charmm type method can be used, but it's
more work than the standard "Amber" type method. It also appears
to be somewhat less efficient in some cases (see Pearlman (1994)
JPC 98, 1487 for more details). In general, I'd stick with the "Amber"
(single topology) method except for cases where you _must_ use
the dual topology (Charmm) approach, e.g. where you have changes
in closed ring topology.
--
David A. Pearlman email: dap@vpharm.com
Vertex Pharmaceuticals Inc.
130 Waverly St. "There are only 25 great people in the
world &
Cambridge, MA 02139-4242 5 of them are hamburgers..." -- Cptn
Beefheart