Subject: MD & MC; Modelling of counterions effects in polyelectrolytes.

From: K Bryson <[email-masked]>
Date: Wed, 6 Jul 1994 19:51:08 +0000

Hi Folks,
		All this recent talk about electrostatic charges on
	protein surfaces and a previous posting about site-specific
	counter-ion binding around the pectic acid in plant cell walls has
	brought up a worry that I've had for some time now, I wonder
	what others think about this.

		This problem applies to all polyelectolyte modelling
	but I will use DNA to explain it since it's what I'm familiar
	with.
		In most of the modelling of DNA we place sodium
	counter-ions next to the -ve charged phosphate groups to
	make the system neutral.
		Now Debye-Manning theory determined that only around
	half the phosphates are neutralised at low ionic concentrations
	and that only when we get to very high salt concentrations do we get
	( I vaguely remember ) around 0.7 of the phosphates neutralised
	by sodium ions, the rest sitting in a counter-ion atmosphere
	some distance away from the DNA.

		So the dilemna in my mind is that we model DNA with
	all the phosphates neutralised very locally as if we are modelling
	it in an extremely concentrated salt solution ( which it would
	denature in ).
		Would it not be better to neutralise half the phosphate
	groups locally with counterions and simulate a charged system ? 
	Statistical mechanics would predict that the DNA would have
	between 0 and 24 ( for dodecamers ) ions bound onto it with
	different probabilities for each depending on salt concentration,
	so why does everyone model it with 24 counter-ions ? When
	at physiological conditions I believe it is most likely to have around
	12 directly bound onto it ?

		This seems to apply to a lot of polyelectrolyte simulations
	 ... full counter-ion neutralisation seems to be what is modelled,
	when electrolytic theory predicts that most polyelectrolytes are
	only partially neutralised by site-specific counter-ions in
	solution.

		What do others think ?
		Maybe I'm being pedantic ( or silly ! ) but the amount
	of counter-ions included in a simulation extensively effects the
	structure of the solute.

		Any better ways to include counterion concentration
	effects ? I suppose a Debye screening term dependant on
	concentration could be applied between all charged atoms separated
	by solvent in the system ... be computational due to the exponential
	term though ... is there a general concensus of belief on this
	problem ?

	Sorry about the large bandwidth.

	Kevin.

=============================================================================
                               
    K.Bryson                 email: [email-masked]   
    Biophysics Group         Tel  : +44 904 430000 Extn. 2236  
    Physics Department       Fax  : +44 904 432214
    University of York      
    Heslington              "Molecular modelling of DNA and its
    YORK, UK                    interaction with small molecules."
    YO1 5DD             

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

---

From: "G. Ravishanker" <[email-masked]>
Date: Wed, 6 Jul 1994 21:11:09 -0400

Kevin

	I would like to comment on the MD simulations of the DNA. We have
spent the past six years trying various methods to simulate the DNA with
and without counterions. As you correctly point out, it is a serious
headache, especially while trying to do a nanosecond level MD on DNA. Here
are some of our observations (We use WESDYN, a derivative of GROMOS, and
the force field used is GROMOS 86 with SPC water. All discussions here
refer to the dodecamer sequence CGCGAATTCGCG):

	1. Our goal is to simulate a "dilute aqueous solution" of DNA.
As a result, we do not apply solute-solute boundary conditions - one
should configure the box/cutoff correctly to avoid some of the artifacts
of arising from this model.

	2. In order for 1 to work correctly, we cannot have the
counterions roam freely, so they are confined to a region around the
surface of the DNA using a harmonic restraint. (In some sense, if we did
this to 75% of the counterions, we will be consistent with the Manning
counterion-condensation model!)

	3. With 1 and 2, the DNA typically looks "OK" (by conventional
convergence criteria) for the first 150ps or so, and then all sorts of
things happen. Basically, the harmonic restraint shifts the ions to be
between two Phosphate groups and the DNA begins to overwind itself etc.

	4. So, we decided to skip the ions altogether. (We, in parallel,
began a NaCl solution simulation to understand the effects of various
protocols on the outcome of the simulation and found really strange
things, such as the like-ions showing an enormous peak in g(R) at whatever
the non-bonded cutoff is. This is likely to happen in the DNA simulation
also. This work was done by Pascal Auffinger and is being currently
written up) We scaled all Phosphate group charges down to -.25 by simply
multiplying the original charges by .25. (CHARMM has a scaled charge
model, where the total Phosphate group charge adds up to -.25, but the
individual atomic charges are still very high)

	5. We have a nanosecond MD with this model, showing very stable
DNA structure. A preliminary communication of this appeared in JACS, vol.
116, 4461-4462 (Authors Kevin McConnell et.al). A detailed paper is under
preparation. 

	We have had reasonable success with this model for various other
sequences ( as measured by the DNA structural stability, and calculated
NMR properties compared against measured properties) also. We are
currently rethinking the whole protocol and we will soon be getting back
to the explicit counterions in some form or the other.

	It also is a fact that for non-canonical B-DNA, the net charges on
phosphates need not be -.24 (see the paragraph below), so anyone wanting
to model Manning's theory, should be constantly readjusting the fractional
charges on the phosphates (may not be much of a change!!!!).

	I offer my sincere apologies, if the following paragraph is a
restatement of a well-known thing and therefore a waste of bandwidth:

	The net charge of -.25 on phosphate groups is valid ONLY for
B-DNA at 298 K. The general condensation formula is:

	f = 1./ (N * Zeta)		Zeta > 1./N
	  = 1.				Zeta < 1./N

	Here f in the uncompensated charge fraction on the phosphate, N is
the valence of the counterion and Zeta is a dimensionless measure of
polyion charge density.

	Zeta = q**2 / (epsilon*k*T*b)

	q is the protonic charge, epsilon, the dielectric constant of
solvent, k , the Boltozmann constant, T, the temperature and b is the
average linear charge spacing between the phosphates along the helical
axis. 

	For an aqueous solution of canonical B-DNA at 298 K, the
fractional charge is -.24 for Na+ and -.12 for Mg2+.

Ravi

****************************************************************************
* Ganesan Ravishanker			Ph: (203) 344-8544 Ext. 3110       *
* Coordinator of Scientific Computing,  Fax:(203) 344-7960                 *
* Adjunct Associate Professor(Dept. of Chem.)                              *
* Wesleyan University               e-mail:[email-masked]  *
* Middletown, CT 06459.                                                    *
****************************************************************************