eef1 (c38b1)
Effective Energy Function 1
EEF1 is an effective energy function combining the CHARMM 19 polar
hydrogen energy function (with certain modifications, see below)
with an excluded volume implicit solvation model. The solvation model
is similar in spirit to the Atomic Solvation Parameter approach, but
does not use surface areas and is therefore much faster. Latest
benchmarks say that simulations with EEF1 take about 50% longer
than the corresponding vacuum simulation.
The solvation model assumes that the solvation free energy of each
group is equal to the solvation free energy of that group in a small
model compound less the amount of solvation it loses due to solvent
exclusion by other atoms of the macromolecule around it. The exclusion
effect of nearest and next-nearest neighbors (1-2 and 1-3 interactions)
are neglected because such neighbors also exist in small model compounds.
The CHARMM nonbonded atom and exlusion lists are used for the solvation
calculation.
Because not only DG but also DH and DCp data are available, we
can calculate the solvation free energy at different temperatures.
This calculation assumes a DCp independent of temperature.
Therefore extrapolation to temperatures very different from 300 K
is not reliable.
EEF1 refers not only to the implicit solvation model but also to
the specific modifications and nonbonded options used in CHARMM.
The nonbonded options must be: ctonnb 7. ctofnb 9. cutnb 10. group rdie
(see example file below).
Three files are needed to use EEF1:
toph19_eef1.inp : This is a modification of toph19.inp where ionic
sidechains and termini are neutralized and contains
an extra parameter type (CR)
param19_eef1.inp: This is a modification of param19.inp which includes
the extra parameter type (CR)
solvpar.inp : This file contains the solvation parameters
When the INTE command is used with EEF1, the number listed under
ASP is the amount of solvation free energy that is excluded between the
two atom selections. For example, the INTE between atom A and atom B will
give the amount of solvation A loses due to B plus the amount B loses due
to A. The command "INTE sele all end" will give the amount of solvation free
energy excluded, not the total solvation free energy of all atoms. That
is, it is not equivalent to "ENERGY".
EEF1 can be used with images. In that case the ASP energy term
refers to the solvation free energy of the primary atoms. This is usually
less negative than when images are not present, because image atoms exclude
some solvation free energy from the primary atoms.
EEF1 is compatible with the BYCC non-bonded option and the NBACtive
command, so that the calculation of non-bonded and solvation energy terms
can be limited to specific subsets of atoms in the system.
The analytical expression of the second derivative matrix of the
EEF1 potential has now been added. Thus, the normal modes, for example,
can now be calculated analytically.
* Syntax | Syntax of the EEF1 commands
* References | Useful references
* Example | Input file
Top
Syntax for EEF1
There are only two EEF1 commands:
EEF1 SETUP [TEMP real] UNIT int NAME solv_param_file
EEF1 PRINT
The first sets up the solvation calculation by giving TEMP
and reading in the solvation parameters. And the second
prints out the solvation of each group. The solvation energy
is stored in ETERM(ASP) and reported under the name "ASP".
Obviously, it makes no sense to use both ASP and EEF1.
If one wants to skip the solvation term after one has set it
up, one can issue the command SKIP ASP.
TEMP is the temperature to which the solvation parameters refer to
(default is 298.15). Note that this is unrelated to the
temperature at which one runs dynamics. It just determines
the solvation free energy parameter values.
PRINT prints out the solvation free energy of each atom/group
as well as the solvation enthalpy and heat capacity
Syntax for EEF1
There are only two EEF1 commands:
EEF1 SETUP [TEMP real] UNIT int NAME solv_param_file
EEF1 PRINT
The first sets up the solvation calculation by giving TEMP
and reading in the solvation parameters. And the second
prints out the solvation of each group. The solvation energy
is stored in ETERM(ASP) and reported under the name "ASP".
Obviously, it makes no sense to use both ASP and EEF1.
If one wants to skip the solvation term after one has set it
up, one can issue the command SKIP ASP.
TEMP is the temperature to which the solvation parameters refer to
(default is 298.15). Note that this is unrelated to the
temperature at which one runs dynamics. It just determines
the solvation free energy parameter values.
PRINT prints out the solvation free energy of each atom/group
as well as the solvation enthalpy and heat capacity
Top
References
[1] T. Lazaridis and M. Karplus, Effective energy function for
proteins in solution, Proteins, 35:133-152 (1999)
[2] T. Lazaridis and M. Karplus, Discrimination of the native from
misfolded protein models with an energy function including
implicit solvation, J. Mol. Biol., 288:477-487 (1999)
[3] T. Lazaridis and M. Karplus, "New View of Protein Folding
reconciled with the Old through Multiple Unfolding Simulations",
Science, 278:1928 (1997)
References
[1] T. Lazaridis and M. Karplus, Effective energy function for
proteins in solution, Proteins, 35:133-152 (1999)
[2] T. Lazaridis and M. Karplus, Discrimination of the native from
misfolded protein models with an energy function including
implicit solvation, J. Mol. Biol., 288:477-487 (1999)
[3] T. Lazaridis and M. Karplus, "New View of Protein Folding
reconciled with the Old through Multiple Unfolding Simulations",
Science, 278:1928 (1997)
Top
---------------------------------------------------------------------
* Example file for EEF1
open read card unit 3 name toph19_eef1.inp
read rtf unit 3 card
close unit 3
open read card unit 3 name param19_eef1.inp
read para unit 3 card
close unit 3
open read unit 3 card name filename.crd
read seque coor unit 3
close unit 3
generate main setup
open read unit 2 card name filename.crd
read coor card unit 2
close unit 2
! IMPLICIT SOLVATION SETUP COMMAND
! The nonbonded options below are part of the model
eef1 setup temp 298.15 unit 93 name solvpar.inp
update ctonnb 7. ctofnb 9. cutnb 10. group rdie
mini abnr nstep 300
!This command prints out solvation free energy for each atom
eef1 print
dynamics verlet timestep 0.002 nstep 1000 nprint 100 iprfrq 100 -
firstt 240 finalt 300 twindh 10.0 ieqfrq 200 ichecw 1 -
iasors 0 iasvel 1 inbfrq 20
inte sele resid 2 end sele resid 19 end
!the command below is not equivalent to energy
inte sele all end
energy
skip asp
energy
stop
---------------------------------------------------------------------
* Example file for EEF1
open read card unit 3 name toph19_eef1.inp
read rtf unit 3 card
close unit 3
open read card unit 3 name param19_eef1.inp
read para unit 3 card
close unit 3
open read unit 3 card name filename.crd
read seque coor unit 3
close unit 3
generate main setup
open read unit 2 card name filename.crd
read coor card unit 2
close unit 2
! IMPLICIT SOLVATION SETUP COMMAND
! The nonbonded options below are part of the model
eef1 setup temp 298.15 unit 93 name solvpar.inp
update ctonnb 7. ctofnb 9. cutnb 10. group rdie
mini abnr nstep 300
!This command prints out solvation free energy for each atom
eef1 print
dynamics verlet timestep 0.002 nstep 1000 nprint 100 iprfrq 100 -
firstt 240 finalt 300 twindh 10.0 ieqfrq 200 ichecw 1 -
iasors 0 iasvel 1 inbfrq 20
inte sele resid 2 end sele resid 19 end
!the command below is not equivalent to energy
inte sele all end
energy
skip asp
energy
stop
Top
New EEF1 parameters (May 2004)
Recent work (e.g. Masunov & Lazaridis, JACS 125:1722,2003) revealed
that the interactions between some ionizable sidechains in EEF1 are too
strong. Also, interactions between hydroxyl groups seem to be too strong.
The files toph19eef1.1.inp and param19eef1.1.inp contain empirical adjustments
of the partial charges to mitigate some of these problems. We refer to this
parameter set as EEF1.1.
In addition, topology files are provided for using EEF1 with the
and top_all22_prot_eef1.1.inp). The standard parameter file can be used with
these. The combination of EEF1 with CHARMM22 has not been extensively tested.
Also, DEBYE-HUCKEL screening of electrostatic interactions has been
implemented, mostly for development purposes. To use it add the keyword
IONIC xxx
where xxx is the ionic strength in mol/lt. With this all electrostatic
interactions are multiplied by exp(-r/rD), where rD is the Debye length
(rD = SQRT(0.0316 Temp/IonicStrength )
New EEF1 parameters (May 2004)
Recent work (e.g. Masunov & Lazaridis, JACS 125:1722,2003) revealed
that the interactions between some ionizable sidechains in EEF1 are too
strong. Also, interactions between hydroxyl groups seem to be too strong.
The files toph19eef1.1.inp and param19eef1.1.inp contain empirical adjustments
of the partial charges to mitigate some of these problems. We refer to this
parameter set as EEF1.1.
In addition, topology files are provided for using EEF1 with the
and top_all22_prot_eef1.1.inp). The standard parameter file can be used with
these. The combination of EEF1 with CHARMM22 has not been extensively tested.
Also, DEBYE-HUCKEL screening of electrostatic interactions has been
implemented, mostly for development purposes. To use it add the keyword
IONIC xxx
where xxx is the ionic strength in mol/lt. With this all electrostatic
interactions are multiplied by exp(-r/rD), where rD is the Debye length
(rD = SQRT(0.0316 Temp/IonicStrength )
Top
Implicit Membrane Model 1
IMM1 is an extension of EEF1 for modeling proteins in lipid membranes
(T. Lazaridis, Proteins, 52:176-92, 2003). The implicit membrane is set up
like this:
open read unit 11 card name toph19_eef1.1.inp
read rtf card unit 11
close unit 11
open read unit 12 card name param19_eef1.1.inp
read para card unit 12
close unit 12
... generate psf, read coordinates ...
eef1 setup membrane slvt water slv2 chex nsmth 10 width 26.0 temp 298.15 -
unit 93 name ../solvpar.new.inp aemp 0.85
... mini, dyna, etc.
The keyword MEMBrane specifies that a membrane is to be modeled. "slvt water"
specifies that the exterior solvent is water and "slv2 chex" that the interior
solvent is cyclohexane. NSMTH (default 10) determines how steep the transition
is at the interface between interior and exterior. WIDTH is the width of the
interior region (default 30A). Standard values are to be used here depending on
the lipid that one wants to model. Such values can be obtained from experimental
data, for example see http://aqueous.labs.brocku.ca/lipid. For example:
DMPC 23.1 A
DOPC 25.4 A
POPC 27.0 A
The last keyword (AEMP, default 0.85) determines the extent of strengthening
of electrostatic interactions in the membrane (the smaller, the stronger).
This parameter was empirically adjusted to give reasonable membrane insertion
energies for model systems.
The above command sets up a neutral/zwitterionic membrane. The effect
of negatively charged lipids can be accounted for by using a Gouy-Chapman term
in the energy function (T.Lazaridis, submitted). This is done by adding the
following keywords:
eef1 setup .... gouy anfr 0.3 area 70. offset 3.0 conc 0.1 valence 1
GOUY specifies that a Couy-Chapman term is to be used. ANFR is the molar
fraction of anionic lipids (e.g., a 70/30 mixture of PC/PG corresponds
to ANFR 0.3, which is the default). AREA is the area (Angstrom^2) per lipid
(default 70). OFFSet is the distance of the plane of negative charge
(usually the phosphates) from the hydrocarbon/water boundary (default 3).
CONC and VALEnce is the molarity and valence of the salt (default 0.1 and 1,
respectively).
CAUTION: When the Gouy-Chapman term is calculated, the ionic sidechains
are given a full charge (they are neutralized otherwise), and this is done
by checking the partial charges. If you want to use topology files other the
ones provided (toph19eef1.1.inp) it might not work.
It is also possible to include the effect of transmembrane voltage
by adding the keyword
VOLT xxx
where xxx is the transmembrane voltage in Volt (default value 0.1).
The transmembrane voltage is set up so that it is positive in the +z
direction. This term is based on the analytical solution to the
Poisson-Boltzmann equation (Roux, Biophys. J, 1997).
The GC and TM voltage energies are added to the Solvation Free Energy
(under ASP column). These terms will be printed out if PRNLEV is greater
than 9.
An implicit cylindrical, toroidal (parabolic) or circular pore in
a neutral membrane can be accounted for by using a modified energy function
(T. Lazaridis, J. Chem. Theory Comput., 1:716-722 (2005); M. Mihajlovic and
T. Lazaridis, Biochim. Biophys. Acta, 1798:1494-1502 (2010)). The following
keywords should be added for the cylindrical, toroidal or circular pore,
respectively:
eef1 setup ... rcyl 10
eef1 setup ... rprb 10 aprb 15
eef1 setup ... rcrc 10
where RCYL specifies the cylindrical pore of 10 Angstrom radius;
RPRB specifies the parabolic pore with the radius in the center of the pore
of 10 Angstrom and APRB defines the curvature of the pore of 15 (which gives
the radius at the ends of the pore of 25 Angstrom); RCRC specifies
the circular pore of 10 Angstrom radius.
Implicit Membrane Model 1
IMM1 is an extension of EEF1 for modeling proteins in lipid membranes
(T. Lazaridis, Proteins, 52:176-92, 2003). The implicit membrane is set up
like this:
open read unit 11 card name toph19_eef1.1.inp
read rtf card unit 11
close unit 11
open read unit 12 card name param19_eef1.1.inp
read para card unit 12
close unit 12
... generate psf, read coordinates ...
eef1 setup membrane slvt water slv2 chex nsmth 10 width 26.0 temp 298.15 -
unit 93 name ../solvpar.new.inp aemp 0.85
... mini, dyna, etc.
The keyword MEMBrane specifies that a membrane is to be modeled. "slvt water"
specifies that the exterior solvent is water and "slv2 chex" that the interior
solvent is cyclohexane. NSMTH (default 10) determines how steep the transition
is at the interface between interior and exterior. WIDTH is the width of the
interior region (default 30A). Standard values are to be used here depending on
the lipid that one wants to model. Such values can be obtained from experimental
data, for example see http://aqueous.labs.brocku.ca/lipid. For example:
DMPC 23.1 A
DOPC 25.4 A
POPC 27.0 A
The last keyword (AEMP, default 0.85) determines the extent of strengthening
of electrostatic interactions in the membrane (the smaller, the stronger).
This parameter was empirically adjusted to give reasonable membrane insertion
energies for model systems.
The above command sets up a neutral/zwitterionic membrane. The effect
of negatively charged lipids can be accounted for by using a Gouy-Chapman term
in the energy function (T.Lazaridis, submitted). This is done by adding the
following keywords:
eef1 setup .... gouy anfr 0.3 area 70. offset 3.0 conc 0.1 valence 1
GOUY specifies that a Couy-Chapman term is to be used. ANFR is the molar
fraction of anionic lipids (e.g., a 70/30 mixture of PC/PG corresponds
to ANFR 0.3, which is the default). AREA is the area (Angstrom^2) per lipid
(default 70). OFFSet is the distance of the plane of negative charge
(usually the phosphates) from the hydrocarbon/water boundary (default 3).
CONC and VALEnce is the molarity and valence of the salt (default 0.1 and 1,
respectively).
CAUTION: When the Gouy-Chapman term is calculated, the ionic sidechains
are given a full charge (they are neutralized otherwise), and this is done
by checking the partial charges. If you want to use topology files other the
ones provided (toph19eef1.1.inp) it might not work.
It is also possible to include the effect of transmembrane voltage
by adding the keyword
VOLT xxx
where xxx is the transmembrane voltage in Volt (default value 0.1).
The transmembrane voltage is set up so that it is positive in the +z
direction. This term is based on the analytical solution to the
Poisson-Boltzmann equation (Roux, Biophys. J, 1997).
The GC and TM voltage energies are added to the Solvation Free Energy
(under ASP column). These terms will be printed out if PRNLEV is greater
than 9.
An implicit cylindrical, toroidal (parabolic) or circular pore in
a neutral membrane can be accounted for by using a modified energy function
(T. Lazaridis, J. Chem. Theory Comput., 1:716-722 (2005); M. Mihajlovic and
T. Lazaridis, Biochim. Biophys. Acta, 1798:1494-1502 (2010)). The following
keywords should be added for the cylindrical, toroidal or circular pore,
respectively:
eef1 setup ... rcyl 10
eef1 setup ... rprb 10 aprb 15
eef1 setup ... rcrc 10
where RCYL specifies the cylindrical pore of 10 Angstrom radius;
RPRB specifies the parabolic pore with the radius in the center of the pore
of 10 Angstrom and APRB defines the curvature of the pore of 15 (which gives
the radius at the ends of the pore of 25 Angstrom); RCRC specifies
the circular pore of 10 Angstrom radius.