CHARMM c32b1 scpism.doc



File: SCPISM ]-[ Node: Top
Up: (commands.doc) -=- Next: Syntax \n


    Screened Coulomb Potentials Implicit Solvent Model (SCPISM)


     The SCPISM is a continuum model of solvated macromolecules that 
treats implicitly the effects of the surrounding aqueous medium. The model
is based on screened Coulomb potentials (SCP) as derived from basic theory 
of polar liquids (the Lorentz-Debye-Sack theory). At the present level
of development the model only incorporates a continuum description
of electrostatic effects and, optionally, a simple term accounting for 
hydrophobic interactions. This implementation is part of a more general 
picture aimed at developing an implicit solvent model (ISM)
for both electrostatic and non-electrostatic effects, to be used in
molecular mechanics simulations of macromolecules. The current 
implementation is suitable for molecular and Langevin dynamics simulations, 
energy evaluation and minimization of peptides and proteins, to be used in 
combination with the all-atom representation (PAR22 or CMAP).
 
     An overview of the physical ideas that led to the development of the
model can be found in http://cmm.cit.nih.gov/~mago/SCPISM.html.
Some technical issues regarding the practical implementation of the SCPISM
into the CHARMM code that are not discussed in previous publications
are also described in the URL above. The present version of the model
will be revised and upgraded in the future. In particular
the parameterization is still subject to development. The parameters used in
all applications of the SCPISM reported to date have not been corrected
or adjusted since the model was first introduced (this is a desirable 
feature of the SCPISM as discussed in the bibliography below, although 
corrections of the present parameters are probably needed, especially those 
governing HB interactions; see [5]). Some ideas towards a parameter-free 
SCP-ISM have recently been reported, which is the ultimate goal for
conceptual and practical necessity [6] (e.g. for polarizable force fields, 
and for transferability of the model to small organic molecules).
In the current implementation there is one parameter per atom type. 
The original parameterization was done based on solvation energies of amino 
Acid side chain analogs and restrictions imposed to the space of parameters,
as reported [2] (the parameterization was not based on reproducing PB
results). Hydrogen bonding strength is treated independently (this was 
shown to be essential for structure prediction and stabilization, as 
discussed in previous publications (e.g. [7])). The refinement of HB 
interactions, along with the introduction of solvent exclusion effects 
and the effect of salt, are currently in development (see [5]).  


* Menu
* Syntax::         Syntax of the SCPISM commands
* Background::     An introduction to the SCPISM (see also URL)
* References::     Useful references (see URL)    
* Example::        Input file 



File: SCPISM ]-[ Node: Syntax
Up: Top -=- Next: Theory -=- Previous: Top \n

    
                          SCPISM commands

An effort was made to minimize the number of input options available 
to the user (i.e., no parameters are allowed to be modified from
outside since the physics of the system is already incorporated into 
the model and, therefore, built within the algorithm). To activate the model 
the following command line is used:

SCPIsm [UISM int]                                                 (1)

where UISM is the unit number for reading the SCP parameters. These
parameters are stored in scpism.inp and must be opened for reading
in the usual way before the model is requested (see example below), i.e.,

OPEN READ UNIT int CARD NAME "scpism.inp"

Once the model is requested, any options for energy, minimization 
and dynamics calculations are supported. Electrostatic 
interactions are truncated using a standard shift function (any other
specification of cutoff of electrostatics is automatically disabled 
once the SCPISM is requested; other options will be incorporated in future 
releases). Only the serial version is supported.
Any restraining option that is introduced as an additional term in the 
potential is supported; use of the CONS FIX constraints are not recommended
because they are inconsistent with the way the self-energies are
calculated: the calculation of Born radii in the SCPISM requires the
positions of all the atoms around each atom of the system to be 
available (in particular the positions of the fixed atoms) because 
it is based on a contact model approach; therefore, it is advised not to
use this CHARMM feature when the SCPISM is on. It is possible (but 
not required) to deactivate the model once any of the required
tasks (energy evaluation, minimization or dynamics) is completed. To 
exit the model use

SCPIsm END

in this case CHARMM will return to the default electrostatics options 
(the model can be switched on again at any point by using command line 
(1) above). Examples for simple MD and LD simulations setup are 
available in http://cmm.cit.nih.gov/~mago/SCPISM.html

Note that, although the current implementation of the SCPISM describes
only electrostatic effects, to be consistent with earlier applications
of the model, an overall hydrophobic term proportional to the total
solvent accessible surface area (SASA), also based on a contact model 
approach, can be requested by using the subcommand HYDRophobic in the
command line (1) above. If this term is not activated, other model
accounting for hydrophobic interactions has to be used.

The CPU time of the SCPISM is shown to be between 2 and 3 times slower
than vacuum (depending on the platform and compilation optimization), for 
cutoff distances between 12.0 and 20 A.

Monte Carlo simulations are not advised with this initial implementation of the 
model because the hydrogen-bonding (HB) algorithm (in particular HB geometry) 
for MD simulations has been simplified (see [2,3,5] and URL above).



File: SCPISM ]-[ Node: Background
Up: Top -=- Next: References -=- Previous: Syntax \n

                  Structure of the parameter input file:

     All the parameters required for the model are atom-type based and are 
collected in a single file. Determination of these parameters was described in 
[2,3].  The parameters have been optimized in the context of the all-atom force 
field.

SCPISM parameter file has the following format (note that not all fields 
in this file are parameters that controls the electrostatics):

H   0.4906  0.6259  0.5000  0.7004  0.3700  0.0052 PH ! polar H
HC  0.5300  9.6746  0.5000  0.7280  0.3700  0.0052 PH ! N-ter H
HA  0.5300  0.5000  0.5000  0.7280  0.3700  0.0052    ! nonpolar H
HT  0.4906  2.3800  0.5000  0.7004  0.3700  0.0052 PH ! TIPS3P WATER HYDROGEN
HP  0.5274  0.5000  0.5000  0.7262  0.3700  0.0052    ! aromatic H
HB  0.4858  0.5000  0.5000  0.6970  0.3700  0.0052    ! backbone H
HR1 0.4651  0.5000  0.5000  0.6820  0.3700  0.0052    ! his he1, (+)his HG,HD2
HR2 0.4580  0.5000  0.5000  0.6768  0.3700  0.0052    ! (+) his HE1

Col. 1: Atom type defined in PAR22
Col. 2: Alpha_i controls slope of D(r) around atom-type i
Col. 3: for PH interaction with PA this value controls the Born radius
        of PH to modulate hydrogen bonding strength (see URL); the
	increase or descrease of the PH Born radius is defined by the
	product Col.3(PH)*Col.3(PA) As a rule, the larger the value of
	this product, the weaker the HB interaction.
Col. 4: Extension of Born radius R_iw to obtain R_ip, i.e.,
        R_ip = R_iw + Col.4 (see [1]); only one value for atoms considered
Col. 5: SQRT(alpha_i); this is needed for alpha_ij=alpha_i alpha_j (see [1,2])
Col. 6: covalent radius for each atom (only 5 values considered, i.e.,
        for C,N,O,S,H)  
Col. 7: gamma_i in hydrophobic energy = summatory over i of gamma_i SASA_i (note
        that all coefficients are equal in this first release) 
Col. 8: denotes atoms involved in H-bonding (PH = polar proton; PA =
        proton acceptor) 


                   Theoretical background of the SCPISM

     The SCPISM is based on screening functions D(r) that modulate the
electric potential phi(r) in the system, rather than the dielectric
functions eps(r) that modulates the electric field E(r). The relation
between both functions is obtained from the definition E(r) = -Grad[phi(r)]

See references and URL above for details.

     Because of the relationship between D(r) and eps(r), both functions
are sigmoidal with r. Once eps(r) is known from theory or experiments,
the screening function that characterizes the medium is obtained by numerical 
integration. Based on these results, the SCPISM uses atom type-dependent 
sigmoidal functions in the context of the all-atom representation PAR22.

     In the SCPISM the standard electrostatic component of the force field
is replaced by terms that describe both the electrostatic interaction energy, 
and the self energy. The screening functions are continuous functions of the 
position and describe a dielectric medium that permeates all of space. For the 
solvated protein, D(r) approaches bulk medium value of the screening only far 
from the protein (see [4] for a discussion). Therefore, the SCPISM does not 
introduce either an internal or an external dielectric constant, and therefore, 
there is no boundary that separates the protein from the solvent.
The SCPISM is derived from basic theory of polar solvation; the model is being 
constructed to incrementally incorporate all the physical effects that are 
removed when the explicit solvent is eliminated. Most important among these 
interactions are non-electrostatic solvent-induced forces (hydrophobic as a 
particular case, but also hydrophilic).

     Because the effective screening functions that characterize the overall
modulation of the electrostatics are obtained from properties of bulk solvent, 
short-range interactions characterizing hydrogen bonding (HB) must be corrected 
(see [3,5] and URL). To obtain a reasonable representation of HB strength in 
solvent medium, all individual HB interactions available in proteins are 
individually calibrated via the self-energy terms (an adjustment of the Born 
radii of polar hydrogens [see [2.3]). The stabilization of HB energies is carried out 
based on charge states of the interacting groups, as well as on hybridization 
states of proton donor and acceptor atoms [5]. For the current implementation 
suitable for MD simulations this approach was dramatically simplified with 
respect to the original development for MC simulations, so that the HB
interactions depend only on atom type (H and HC interacting with O, OC, OH1
and NR2) and no directionality has been included (this was shown to be important,
however, for structure calculation based on MC simulations to avoid artifacts
such as multiple HB between donor and acceptor) 



File: SCPISM ]-[ Node: References
Up: Top -=- Next: Example -=- Previous: Background \n


                               References


[1] S A Hassan, E L Mehler, D Zhang and H Weinstein, Molecular Dynamics
    Simulations of Peptides and Proteins with a Continuum Electrostatic Model
    based on Screened Coulomb Potentials; Proteins 51, 109 (2003)

[2] S A Hassan, F Guarnieri and E L Mehler, A General Treatment of Solvent 
    Effects based on Screened Coulomb Potentials; J Phys Chem B 104, 6478 
    (2000)

[3] S A Hassan, F Guarnieri and E L Mehler, Characterization of Hydrogen 
    Bonding in a Continuum Solvent Model; J Phys Chem B 104, 6490 (2000)

[4] S A Hassan and E L Mehler, Critical Analysis of Continuum Electrostatics:
    the Screened Coulomb Potential Implicit Solvent Model and the Study of
    the Alanine Dipeptide and Discrimination of Misfolded Structures of
    Proteins; Proteins 47, 45 (2002)

[5] S A Hassan, Intermolecular Potentials of Mean Force of Amino Acid Side Chain    
    Interactions in Aqueous Medium; J Phys Chem B 50, 19501 (2004)

[6] S A Hassan and E L Mehler, From Quantum Chemistry and the Classical
    Theory of Polar Liquids to Continuum Approximations in Molecular Mechanics
    Calculations; Int. J. Quant. Chem. 102, 986 (2005)

[7] X Li, S A Hassan and E L Mehler, Long Dynamics Simulations of Proteins using 
    Atomistic Forcefields and Continuum Representation of Solvent Effects: 
    Calculation of Structural and Dynamics Properties; Proteins X, XX (2005) 


File: SCPISM ]-[ Node: Example
Up: Top -=- Next: Top -=- Previous: References \n

                 
                         SCPISM example input file

---------------------------------------------------------------------
* energy, minimization and dynamics with the SCPISM
*

open read card unit 10 name top22.inp
read rtf unit 10 card
close unit 10 

open read card unit 10 name par22.inp
read para unit 10 card 
close unit 10 

open read unit 10 card name SCPISM.inp

! define system
generate main setup

open read unit 2 card name filename.crd
read coor card unit 2 
close unit 2

! set up all options for the run (nbond cutoff distance, shake, vdW, etc)
! these options can also be specified in the 'energy' command line, as usual

! request SCPISM 

SCPI HYDR UISM 10

! this will activate electrostatics and a simple hydrophobic term;
! if other model is used to describe hydrophobic interactions then
! replace the above command line by: SCPI UISM 3

! now calculate energy, minimize structure and run dynamics

ener 

mini [options] 

dyna [options]

! exit SCPISM

SCPI END

! if needed, the model can be requested again at any point

stop
-----------------------------------------------------------------------


CHARMM .doc Homepage


Information and HTML Formatting Courtesy of:

NIH/DCRT/Laboratory for Structural Biology
FDA/CBER/OVRR Biophysics Laboratory
Modified, updated and generalized by C.L. Brooks, III
The Scripps Research Institute