Article Outline

Article Information

PubMed

Related Articles

  • …more

Copyright © 1999 The Biophysical Society. All rights reserved.
Biophysical Journal, Volume 77, Issue 4, 1858-1870, 1 October 1999

doi:10.1016/S0006-3495(99)77029-3

Biophysical Theory and Modeling

Brownian Dynamics Simulation of DNA Condensation

Pierre-Edouard Sottas, Eric Larquet, Andrzej StasiakGo To Corresponding Author,   and Jacques Dubochet

Laboratoire d’Analyse Ultrastructurale, Bâtiment de Biologie, Université de Lausanne, CH-1015 Lausanne, Switzerland

Address reprint requests to Dr. Andrzej Stasiak, Lab. d’Analyse Ultrastructurale, Bâtiment de Biologie, Université de Lausanne, Niveau 1, CH-1015 Lausanne-Dorigny, Switzerland. Tel.: 41-21-692-4282; Fax: 41-21-692-4105.

Abstract

DNA condensation observed in vitro with the addition of polyvalent counterions is due to intermolecular attractive forces. We introduce a quantitative model of these forces in a Brownian dynamics simulation in addition to a standard mean-field Poisson-Boltzmann repulsion. The comparison of a theoretical value of the effective diameter calculated from the second virial coefficient in cylindrical geometry with some experimental results allows a quantitative evaluation of the one-parameter attractive potential. We show afterward that with a sufficient concentration of divalent salt (typically ∼20mM MgCl2), supercoiled DNA adopts a collapsed form where opposing segments of interwound regions present zones of lateral contact. However, under the same conditions the same plasmid without torsional stress does not collapse. The condensed molecules present coexisting open and collapsed plectonemic regions. Furthermore, simulations show that circular DNA in 50% methanol solutions with 20mM MgCl2 aggregates without the requirement of torsional energy. This confirms known experimental results. Finally, a simulated DNA molecule confined in a box of variable size also presents some local collapsed zones in 20mM MgCl2 above a critical concentration of the DNA. Conformational entropy reduction obtained either by supercoiling or by confinement seems thus to play a crucial role in all forms of condensation of DNA.

Introduction

Condensation of linear DNA in aqueous solution requires a sufficient concentration of multivalent counterions. Since within viruses and cells DNA generally adopts a highly compact and condensed structure requiring charge neutralization of DNA, the action of such multivalent cations has been intensively studied (see for review Bloomfield, 1996,Bloomfield, 1997). DNA condensation studies are also motivated by biomedical applications to gene therapy since the discovery that introduction of DNA into host cells is facilitated by the presence DNA precipitates formed in the presence of divalent salt (Felgner and Ringold, 1989).

Condensation of linear DNA in vitro requires counterions with a valence of at least +3 in an aqueous environment (Bloomfield, 1996). In vivo solutions contains a large amount of univalent and divalent salts like KCl, NaCl, and MgCl2, and also organic polyvalent counterions like spermine or spermidine. The level of DNA condensation in living cells is regulated by changes of the ratio between concentrations of mono, bi, and polyvalent counterions. However, even at constant ionic conditions DNA condensation can be affected by changes of a torsional stress in DNA.

A few years ago, Bednar et al and Ma and Bloomfield, 1994 provided evidence that concentrations of divalent counterions near the physiological ones can lead to a condensed form of supercoiled DNA, whereby torsionally relaxed DNA does not condense under the same conditions. Supercoiling appears therefore to cooperate with the screening effect of the counterions in enhancing the close lateral association of DNA helices required for condensation. These observations show that torsional stress can promote DNA condensation by lowering the conformational entropy of DNA molecules.

A Brownian dynamics (BD) simulation of polymers (Klenin et al) is a suitable way to investigate DNA dynamics, and we decided to use this method to study the behavior of the plasmid pUC18 (2686 bp) upon varying concentrations of counterions, level of supercoiling, and concentration of DNA. This work is also motivated by the great progress in recent years in polyelectrolyte theory, where it is now demonstrated that divalent counterion fluctuations lead to the emergence of an attractive electrostatic potential (Bhuiyan and Outhwaite, 1994,Gronbech-Jensen et al,Gronbech-Jensen et al).

Counterion-mediated DNA condensation cannot be explained by just a screening effect of the phosphate charges present on the DNA molecule. This effect can only reduce the distance between opposite segments of the double helix and lead to a tighter supercoiled form of the circular DNA without a transition to a condensed form, where opposing segments of interwound superhelices show zones of lateral contact with each other. Without condensation, the superhelix diameter (Vologodskii and Cozzarelli, 1994) decreases continuously in increasing salt solutions and always remains much larger than the geometrical diameter of B-DNA. The effect of screening by counterions is a pure electrostatic effect and is well described at low salt by the mean-field Poisson-Boltzmann (P-B) approximation (Israelachvili, 1992,Hecht et al). This formalism is integrated in the simulation in a form first suggested by Stigter (1977; Stigter and Dill, 1993) and extended here to divalent salt solutions. As two parallel cylinders with the same line charge always should repel each other in this approximation, the explanation of the condensation of polyelectrolytes needs the addition of an attractive potential. The condensation mechanism is not thus primarily electrostatic.

The origin of attractive forces between macromolecular bodies mediated by a solvent is an old and well-known problem. The double-layer theory in colloidal physics based on the DLVO theory (Israelachvili, 1992) assumes a direct additivity of a repulsive potential and an attractive one, originating in van der Waals forces, which has in many cases underestimated the size of the attractive part. In the last few years some new theories and simulations at a molecular level (Gronbech-Jensen et al,Gronbech-Jensen et al,Das et al,Gavryushov and Zielenkiewicz, 1998) lead to a revival of the question of the nature and the features of the attractive potential. Attractive free energy can take its origin in London-like dispersion forces between opposite segments of the DNA molecule (Bloomfield et al,Israelachvili, 1992), in the counterion fluctuations (Oosawa, 1968,Bhuiyan and Outhwaite, 1994,Rouzina and Bloomfield, 1996) and in the reorganization of the water molecules surrounding the DNA surface (Rau and Parsegian, 1992).

Because all these attractive forces act on length scales very small and similar to each other, it is very difficult to disentangle them. Furthermore, the BD of polymers with hydrodynamic interactions presented here deals only on a mesoscopic scale with a characteristic length on the order of several Debye lengths, and does not allow specification of the nature of the attractive interactions that act on atomic dimensions. Because at this scale the studied system is very complex due to the high number of degrees of freedom and interactions, we think that only an approach on an atomic level, such as molecular simulations (Sprous and Harvey, 1996,Gavryushov and Zielenkiewicz, 1998), will give a correct explanation of the nature of these forces. Despite the fact that much progress has been made theoretically, no theory can nowadays model these forces in the geometry of two interacting double-helices of DNA. Recent developments in polyelectrolyte theory with the introduction of the modified Poisson-Boltzmann (MPB) equation, and taking into account both the charges fluctuations and the ionic exclusion volume effects (Bhuiyan and Outhwaite, 1994), have been used for the calculation of the electric potential around DNA (Gavryushov and Zielenkiewicz, 1998). These last studies predicted an attractive potential around DNA even for diluted divalent electrolytes. The MPB equation can thus be applied to studies of the electrostatic contribution to DNA-DNA interactions in the presence of multivalent ions or at relatively high concentrations of monovalent ions.

In this work we address the question of a qualitative interpretation of the attractive potential and try to quantify this potential by determining the second effective virial coefficient of DNA in divalent salt solutions. The overall attractive effects are accounted for by a simple van der Waals potential between cylinders with only one adjustable parameter distinguishable as an effective Hamaker constant.


Theory

The DNA condensation arises from a complex interplay of interactions (Bloomfield, 1997). In order to introduce all these interactions into a manageable simulation, we need to assume that each of them is disentangled from the other and is determined by only one parameter. Our simulations of BD take into account the wormlike chain model of the polymer, the electrostatic repulsions screened by the counterions, the effect ascribed to attractive forces, and the repulsive excluded volume due to the steric properties of DNA. Finally, the theory of the second virial coefficient of DNA is presented, which takes into account these three types of potentials in a cylindrical symmetry and allows a quantitative evaluation of the one-parameter attractive potential.

Brownian dynamics of polymers

To express the global features of the structure of DNA, including time-dependent transitions, a classical approach by Monte Carlo simulation is not sufficient and has to be replaced by a BD simulation with hydrodynamic interactions. This approach of describing the dynamics of a polymer includes the thermal motion of the surrounding molecules of the solvent and has already been applied with success to linear and circular DNA (Allison et al,Chirico and Langowski, 1994,Klenin et al,Jian et al,Merlitz et al).

The basic idea of the BD algorithm is to present the DNA molecule as a bead-chain where the hydrodynamic interactions with the assumed homogeneous viscous solvent are treated using the Rotne-Prager tensor. The expression of the thermal motion of the system is done via the introduction of a random Langevin force in the equations of motion with a discrete characteristic time of typically 50ps. The integration of these equations for the three degrees of freedom of beads position (Ermark and McCammon, 1978) and for the fourth degree of the torsional motion of the polymer gives the basic algorithm of the kinetics of the system. Klenin et al have recently proposed a general-purpose BD program for computing polymer dynamics. In this work we have followed, in large part, their model except for the computation of the electrostatic interactions that are presented separately in the next section. In brief, the DNA molecule is modeled as a wormlike polyelectrolyte chain with both electrostatic and hydrodynamic interactions. The chain consists of N straight segments that are connected through appropriate bending, torsion, and stretching forces derived from harmonic potentials. The exact expressions for these mechanical forces can be found in Klenin et al. The bending elasticity of the polymer in the absence of electrostatic interactions is derived from the DNA persistence length known to be of ∼50nm. Electrostatic interactions are known to affect the effective persistence length in a range of 10% depending on salt conditions (Vologodskii and Cozzarelli, 1995), and this also occurs in our simulations due to the contributions of the electrostatic forces. In all the computations referred to in this paper, the torsional persistence length is set to 63nm (corresponding to 2.6 · 10−19 erg · cm at ambient temperature), which is a compromise between several experimental values ranging from 50 to 75nm. The stretching stiffness is chosen to have a small variance of the average segments length lo (lo2/300). The value of this parameter does not influence the dynamics significantly, as already pointed by Klenin et al, but smaller values of the stretching stiffness allow larger timesteps (Jian et al). The typical molecule studied in this work is the pUC18 plasmid (2686 bp) and is represented by 178 segments, each having a mean length of 5nm. This gives 20 segments per Kuhn statistical length (or 10 segments per persistence length), a value known to guarantee good statistical properties of the simulated chains (Vologodskii and Cozzarelli, 1995). Because the segment's length is small enough compared with the persistence length, the corresponding DNA segment can be considered rigid in a first approximation and viewed as a cylinder of radius equivalent to the geometrical radius of the B-DNA molecule (1nm). However, in order to model the hydrodynamic interactions using the Rotne-Prager tensor as proposed in Ermark and McCammon, 1978, we have to abandon this cylindrical symmetry of the wormlike chain model in attaching a spherical bead at each vertex of the chain. The connection between these two approaches is possible via the definition of the bead radius. This radius is calculated as in Klenin et al in such a way that the hydrodynamic properties of a chain of beads are the same as the hydrodynamic properties of a polymer of the same length. This is done in equaling the diffusion coefficient of these two entities with a hydrodynamic radius of the cylinder equal to 1.2nm. This radius is a little larger than the geometrical radius of B-DNA to account for the first hydration shell composed of fixed molecules. The radius of the phantom beads obtained by this way then allows the computation of the polymer dynamics according to the Ermark-McCammon algorithm. Several other environmental parameters, such as the temperature (293K), the dielectric constant of water (80), and the solvent viscosity (1 Cp) are kept constant during all the simulations. The dynamics timestep was chosen to have a typical single displacement of a bead on the order of several angstroms. This length is small enough to explore the phase space of two interacting DNA cylinders. Thus a value of 50ps was set when attraction interactions were small and decreased to 25ps in the presence of methanol. Larger values of timestep are avoided because they can lead, in case of strong attractive interactions, to a change of the topology with a crossing of two segments during the simulation.

The evaluation of the hydrodynamic interactions is the most time-consuming part in this BD simulation since it includes a Cholesky factorization of a (3N) · (3N) matrix, and the small value of the timestep chosen here can discourage us from studying long dynamics. However, the simulation can be accelerated significantly by using the same hydrodynamic matrix over several timesteps instead of recalculating it step by step. Klenin et al have shown that a 100-step update interval does not significantly change the dynamical properties of circular and linear noncollapsed DNA. As we can think that the dynamics of two lateral segments with a separation on the order of 3nm is somewhat more sensitive, here we use only a 10-step update of this matrix to be on the safe side.


Electrostatic potential of DNA in salt solutions

The parametrization of the immediate ionic environment of a polyelectrolyte molecule like DNA needs two assumptions. First, because a numerical computation at the atomic level in a mesoscopic system is unthinkable nowadays, the DNA is idealized as a uniformly charged cylinder. Close to DNA the cylindrical symmetry is broken, but this approximation remains relevant at distances on the order of several Debye lengths. The second assumption deals with the integration of the P-B equation. This nonlinear differential equation has recently found an analytical solution for objects with cylindrical symmetry (Tracy and Widom, 1997), but this exact solution is unfortunately very time-consuming for this type of computer simulation. However, as the potential that prevails near the polyelectrolyte molecules is very high due to the superposition of the fields of a large number of charges, the Debye-Hückel (D-H) approximation, which leads to a linear form of the P-B equation, is irrelevant. In this work we follow an intermediate approach originally proposed by Stigter, 1977.

The double layer is divided into three regions. From a region inside the surface of hydrodynamic shear of DNA called the Stern layer starts a domain where high potentials take place and the nonlinear equation prevails. Finally, the high potential curve approaches a D–H curve at large distance from the cylinder. The major advantage of this model is that the linearized P-B equation is solvable in closed form in the third region and can be extended numerically into the nonlinear domain. The first region is not delimited by the steric radius of DNA, but rather by a more realistic cylindrical surface that separates mobile ions on the outside from intrinsic DNA charges plus fixed counterion charges on the inside. Thus, as for the computation of the hydrodynamic interactions, this cylindrical surface is chosen with a radius xo=1.2nm.

The P-B relation between the potential φ and the charge density ρ is

(1)
(2)
where D is the dielectric constant of the solvent, Ni the number of elementary charges q of species i, and Ci the concentration of ions i far from the polyelectrolyte molecule.

Introducing the Debye length 1/κ as the unit of length and a dimensionless potential ψ given by

(3)
one gets for 2:1 salt as MgCl2 the P-B equation in cylindrical coordinates
(4)
For ψ≪1 that is for large x, Eq. (4) can be linearized and the D-H equation is obtained
(5)
A solution of Eq. (5) is given by the zeroth-order modified Bessel function of the second kind, K0(x). The integration over the distance of this solution leads to the P-B potential between two straight thin lines
(6)
where R is their shortest distance and θ their angle of alignment.

The basic idea of Stigter is to use an effective rather than actual charge density of the polyelectrolyte (Stigter, 1977). Thus, the solution of the D-H equation for a straight thin line with effective charge density γeff should coincide in the third region with the solution of the P-B equation for a cylinder of charge density γ.

Another problem is a suitable definition of the Stern layer. As the atomic dimensions prevail in this domain, the understanding of its structure and energetics requires the evaluation of the interactions of ions of finite size as regulated by the topography of the DNA-water interface. The cylinder is modeled as uniformly charged and, to take into account the immobile counterions inside the shear surface, the intrinsic charge density of the polyelectrolyte is thus lowered by a factor α. This procedure is known to give a realistic representation of the electrostatic potential at low salt (Stigter and Dill, 1993; Schlick et al,Hecht et al) and spares computer time.

For 1:1 salt solutions, the counterions have a low affinity to bind directly with the phosphate group and are in a state of essentially complete hydration and free translational and rotational mobility. However, it is known that divalent cations interact strongly with DNA (Duguid et al). Consequently, the factor α expresses the quasi-fixed charges inside the Stern layer. Over the last two decades a great deal of work was done on DNA-metallic cations interactions, and α values can be found directly with the calculation of the adsorption isotherms of counterions in the Stern layer. Experimental results can be found in Skerjanc and Strauss, 1968, Clement et al, Record et al, and Li et al. Note that for 1:1 salt solutions α remains almost constant in a large range of concentrations, while for 2:1 salt solutions α decreases with increasing concentrations of salts. This indicates that some water molecules are released from the first hydration shell and are replaced by some divalent counterions.

One gets at last the linear charge of the modeled cylinder

(7)
and applying Gauss's law in cylindrical coordinates, we have
(8)
The effective charge density γeff can now be found numerically by integrating Eq. (4) with boundary conditions given by Eqs. (5).

Finally, the electrostatic energy between two cylinders of DNA in the simulation is approximated by a pairwise summation of linear D-H screening potentials for single charges along the segments li

(9)
where lb is the Bjerrum length and Uij is given in units of kbT.

During the simulation procedure, Eq. (9) is approximated as in Schlick et al and Vologodskii and Cozzarelli, 1995 with a subdivision of the segment in m parts of equal distance and the derivatives of the potential are calculated numerically. An exact integration of Eq. (9) in the BD algorithm is given in Klenin et al but needs a preliminary computation of the derivatives of the potential. Table 1 gives the results of this approach for a large range of concentrations of NaCl and MgCl2.

Table 1 Electrostatic data for DNA cylinders in the mean-field P-B approximation
Concentration (M)Inverse Debye length κ (nm−1)Effective diameter db (nm)Charge correction factor αHydrodynamic radius x0 (nm)(nm−1)
NaCl
13.282.970.731.24.10E+03
0.52.323.420.731.2482
0.21.474.410.731.271.5
0.11.045.620.731.226.5
0.050.7367.420.731.212.8
0.020.46511.20.731.26.27
0.010.32915.70.731.24.24
0.0050.23322.30.731.23.13
MgCl2
0.11.82.220.081.22.27
0.051.272.730.141.22.33
0.020.8044.140.231.21.60
0.010.5685.950.301.21.22
0.0050.4028.550.371.20.917
0.0020.25413.90.471.20.646
0.0010.1820.30.521.20.505

Attractive potential of DNA

Condensation of linear DNA with increasing salt concentration of valence +3 is a well-known fact; however the mean-field P-B theory, which predicts a purely repulsive force at all separations for all valences, fails to explain this transition. Here one assumes that even at high concentration of salt of all valences such a theory gives a good representation of the repulsive interactions and that condensation can only be explained by the addition of an attractive potential.

In the last few years, several groups have attempted to describe the attraction between highly charged macroions in multivalent salt solutions (Rau and Parsegian, 1992,Bhuiyan and Outhwaite, 1994,Rouzina and Bloomfield, 1996). The exact origin of these forces is not yet well-known, and only an atomic level approach (Gronbech-Jensen et al,Gronbech-Jensen et al,Gavryushov and Zielenkiewicz, 1998) will lead to a better understanding of the problem.

It is well-known that ordinary van der Waals forces between segments of a double helix cannot induce DNA condensation by themselves. In a quantitative way, counterion fluctuations, hydration forces, and probably in a weaker part, discreteness of DNA charges contribute to attractive effects. However, an attempt to disentangle them presents obvious difficulties. That is why we model the attractive potential as coming from simple effective van der Waals forces between cylinders of macromolecular dimensions.

In the early 1970s, several authors attempted to calculate the van der Waals energy between cylinders (Parsegian, 1972,Brenner and McQuarrie, 1973). Here, we approximate the attractive potential between two cylinders of infinite length and inclined by an angle θ by:

(10)
where the Di are numerical constants (Brenner and McQuarrie, 1973), r is the hydrodynamic radius of the DNA molecule (1.2nm), R the shortest distance between their axes, and A the adjustable parameter. A can be seen as an effective energy and has no direct relation with the Hamaker constant calculated from classical Lifshitz theory on dispersion forces.

During the simulation, the attractive energy is computed as follows. The energy between two nonadjacent cylinders is approximated by Eq. (10) with a cutoff at 10nm, beyond which it is assumed that the attractive part becomes negligible. To avoid numerical divergence for small angles between the rods, the factor 1/sin(θ) is replaced by a smaller number of the order L/R, where L denotes the length of the cylinder. That is when θ<0.1 we use

(11)
which tends asymptotically for large R to the well-known van der Waals potential between parallel cylinders with a dependence in 1/r5.


Theory for the second virial coefficient of DNA

Another way to express the electrostatic interactions is in terms of DNA effective diameter db (Onsager, 1949). This parameter is the only characteristic of an idealized hard-core potential that mimics the impenetrability between polymers. A nice feature of this parameter is its simplicity and that its theoretical value, given by the second virial coefficient, can be directly compared with various experimental methods. The second virial coefficient describes the thermodynamic properties of the pair interaction between solute molecules. In the case of two interacting cylinders, the second virial coefficient B2 is given by:

(12)
where U(r1, r2) is the pair interaction potential between the two rods in position r1 and r2 in the volume V. In case of a hard-core potential parametrized by db, it becomes
(13)
This result allows us to define an effective diameter from a continuous theory of polyelectrolytes. For charged rods in monovalent solutions, db is calculated by introducing Eq. (6) into Eq. (12), and is well-approximated by
(14)
where F depends on the effective charge density of DNA γeff (Onsager, 1949). For a high concentration of divalent salt, the relation in Eq. (14) does not hold anymore due to the lower effective charge density of DNA, and a numerical integration of Eq. (12) is necessary.

The values of the DNA effective diameter calculated from the P-B equation in NaCl solutions (see Table 1) can be compared with various experimental data (Vologodskii and Cozzarelli, 1995) and give an excellent agreement. For MgCl2 salt, fewer experimental data are available but comparison with studies of DNA knots formation by Shaw and Wang, 1993 and Rybenkov et al is possible. As shown in Fig. 1 a very good agreement between theoretical and experimental values is obtained at low salt concentrations (2mM and 5mM MgCl2). It's not surprising to find such a good correspondence at these concentrations due to the mean-field approximation of the P-B theory. However, at higher concentrations, the P-B theory always predicts repulsion, although the experience gives an effective diameter lower than the geometric diameter of B-DNA starting from ∼20mM MgCl2, with a discrepancy increasing with the concentration (Shaw and Wang, 1993,Rybenkov et al). At this range of concentrations, attractive forces, especially forces arising from the counterion fluctuations, cannot be neglected.

Display large version of this figure
Figure 1
Effective diameter db of DNA as a function of MgCl2 concentrations. Data along dashed lines are from Shaw and Wang, 1993 and were obtained by measuring knotting probability upon circularization of 5.6kb (○) and 8.6kb (△) linear DNAs. The data were adjusted with the formula proposed by Rybenkov et al. The horizontal solid line corresponds to the geometric diameter of B-DNA. The isolated experimental point (▿) is from Rybenkov et al at 20mM MgCl2 and 1mM NaCl. The solid line is calculated from the renormalized Poisson-Boltzmann mean-field theory.

Another approximation arises now from the fact that experimental processes cannot determine the third virial coefficient, whereas counterion fluctuations are typical features that fall into this coefficient (Post and Zimm, 1982). Thus the second virial coefficient obtained experimentally is seen as an effective one and the attractive potential is directly introduced in its definition. Accepting this assumption, we will now focus on the 20mM MgCl2 concentration, but it is straightforward to extend what follows to other concentrations. There are several reasons to fix our attention on this concentration. First of all, this is exactly the concentration at which the effective diameter equals the steric diameter of DNA. Second, it is known that at this concentration linear DNA does not condense in aqueous solutions, whereas the introduction of 50% methanol produces an aggregation. Third, we know that polymers present very interesting properties exactly at the Flory θ-point where the second virial coefficient is also characterized by a cancellation of the repulsive electrostatic forces by attractive ones. Finally, this concentration gives an ionic force close to the physiological one in a high number of biological specimens.

The introduction of the attractive potential in the second effective virial coefficient allows us to quantify the attractive parameter A introduced in Eq. (10). This effective Hamaker constant is thus found by equaling the effective diameter obtained from experimental data and the integration of Eq. (12) while accounting for the attractive potential given by Eq. (10). The impenetrability of two DNA cylinders is parametrized by an exp(−r/λ) potential, which is typical of a repulsive hydration potential, with λ the decay rate (0.3nm, Rau and Parsegian, 1992). The results are not too sensitive to the form of this steric potential but rather to the distance from which it occurs. We thus introduce a cutoff of the potential at 2.6nm. This length is approximately the interhelical spacing in condensed DNA induced by Mg2+ in 50% methanol solutions obtained from x-ray diffraction studies (Schellman and Parthasarathy, 1984).

Note that for 20mM MgCl2 the effective diameter is very close to the geometric diameter of B-DNA (see Fig. 1), that is, the effect of the attractive potential makes up exactly for the effect of the repulsive part in the calculation of the second virial coefficient. Furthermore, both P-B and attractive potentials have the same angle's dependence in 1/sin(θ) (Eqs. (6)) and it is sufficient for fixing the attractive parameter to have the relation

(15)
where W(R) is the radial part of the total interaction potential. The effective attractive energy parameter takes in this formalism the value A=0.86 kbT in 20mM MgCl2. Bloomfield et al have described a similar method whereby a value between 2 kbT and 5 kbT can produce a net attraction at a separation of ∼3nm at a similar ionic force. The discrepancy between the two studies comes from the fact that we have chosen a hydrodynamic radius equal to 1.2nm instead of 1nm. A check with this smaller radius gives, in our formalism, an effective Hamaker constant of 2.1 kbT, which is in agreement with Bloomfield et al previous calculations.

The potential W(R) is plotted in Fig. 2 in comparison with the potential obtained from the mean-field P-B theory. W(R) can be seen as the energy between two perpendicularly crossed cylinders of DNA of infinite length. It is very interesting to compare the values of the maxima of the potential obtained from this theory and experimental data from osmotic pressure measurements (Rau and Parsegian, 1992). These measurements did not indicate a clear transition toward a condensed form of DNA in the presence of 20mM MgCl2 and 10mM TrisCl. However, these studies of DNA interactions indicated that the DNA interaxial spacing is almost independent of the condensing agent and a transition can occur from ∼3.2 to ∼2.8nm (Rau and Parsegian, 1992). These values can be compared with the positions of the maxima of the radial potential W(R), respectively, 3.1 and 2.6nm.

Display large version of this figure
Figure 2
Radial electrostatic potential between DNA helices at 20mM MgCl2. The dashed line corresponds to the renormalized Poisson-Boltzmann theory. The solid line represents the electrostatic potential with the incorporation of a short-range attractive potential and a repulsive one accounting for the steric properties of DNA. This potential gives an effective diameter equal to the geometric diameter of B-DNA (2nm).


Results

The DNA condensation threshold depends on the force of the attractive potential and on the difference of the entropy of the system before and after the condensation. DNA supercoiling strongly reduces the conformational entropy of DNA molecules and thus can favor DNA condensation. As discussed in the Introduction, earlier experimental studies demonstrated that ∼20mM MgCl2 is sufficient to bring down the effective diameter of DNA to the same value as its geometrical diameter. We therefore decided to simulate the effect of 20mM MgCl2 on supercoiled DNA molecules having a supercoiling density of −0.05, which is close to a natural one in DNA isolated from Escherichia coli cells. To demonstrate the importance of attractive potential contributions to the mean-field P-B electrostatic potential we compare simulations neglecting them with these accounting for them. Then we show that when attractive potential is very strong, as in 50% methanol, one observes DNA condensation without requiring DNA supercoiling. As a negative control of DNA condensation, we show that in aqueous solutions with 20mM MgCl2 the DNA does not condense when it is not supercoiled. Finally, we show that other ways of reducing conformational entropy, such as confinement to a reduced volume, can also favor DNA condensation similarly to DNA supercoiling.

The program was implemented in Fortran90 and executed on two DecAlpha Workstations, 400au and 500au, respectively; to check the algorithm and to have sufficiently accurate statistics, >5 months of CPU time on these two machines was needed.

Supercoiled DNA in 20mM MgCl2 in simulations neglecting attractive interactions

To characterize the shape of supercoiled DNA molecules a good measurable parameter is the superhelix diameter, which tells us what the local distance is between two opposing segments of the supercoiled DNA molecule. The distribution of this diameter at 20mM MgCl2 obtained in simulations neglecting the attractive electrostatic potential is shown in Fig. 3. The size of the DNA (2686 bp) corresponds to the pUC18 plasmid. The supercoiling density is chosen near −0.05 (ΔLk=−13). The other parameters of the simulation can be found in the section devoted to the description of the BD algorithm. This distribution is obtained after an equilibration time of ∼0.2ms and measured for several randomly chosen configurations. The inset of Fig. 3 presents the evolution of the writhe in these conditions and demonstrates that 0.2ms suffices to reach equilibrium. It is an accepted fact that thermal ensembles generated by a long BD run correspond to the equilibrium distribution of conformations. Note that the most frequent distance between opposing segments of plectonemically interwound DNA molecules is close to 7nm, and that all the values are larger than the geometrical diameter of DNA. A similar slope was found by Hammermann et al by Monte Carlo simulations of DNA in elevated NaCl concentrations under conditions where repulsion is screened but where there is no attraction. High values of the superhelix diameter measurable over certain positions of simulated molecules are typical features of the endloops. Simulated molecules are frequently unbranched or contain one branching point of their superhelical axis (a few cases with more branching points were also seen, data not shown). A relatively low ratio of branched molecules can be explained by a lack of intrinsic curvature in modeled DNA molecules. Sprous and Harvey, 1996 demonstrated, using molecular dynamics simulations, that branching is promoted by the presence of intrinsically curved regions. The number of branches does not seem to affect the superhelix diameter measured outside endloops. A typical conformation obtained at equilibrium and having one branching point is shown in Fig. 4.

Display large version of this figure
Figure 3
Distribution of the superhelix diameter at 20mM MgCl2 without the introduction of an attractive potential in the simulation. The superhelix diameter was calculated as the minimal length between a subunit and its nearest segment while excluding a suitable number of nearest neighbors. Loops are note excluded in this calculation. The distribution was calculated over several independent conformations of a DNA of 2686 bp with a supercoiling density equal to −0.05. Each configuration was obtained after a preliminary run of 0.2ms from when it is assumed that the equilibrium is obtained as shown by the evolution of the writhe in the inset.
Display large version of this figure
Figure 4
Supercoiled DNA of 2686 bp (pUC18) with a difference of linking number equal to −13 in 20mM MgCl2 without the introduction of an attractive potential into the simulation.

Supercoiled DNA in 20mM MgCl2 in simulations taking into account the attractive interactions

The introduction of the attractive potential as described in the Introduction strongly changes the shape of the molecule. The obtained conformations now present some zones of lateral contacts between segments, as we can see in Fig. 5 even at 20mM MgCl2 in aqueous solvent. All the parameters of the simulation were kept the same as before except for the attractive potential parameter, which now has a non-zero value (A=0.86 kT). The corresponding distribution of the superhelix diameter over several configurations is shown in Fig. 6. The maximum at ∼2.5nm characterizes the collapsed zones. However, there is a surprising fact that collapsed zones are frequently interspersed with open zones. In defining a condensed zone as a zone where the superhelical diameter is below 4nm and integrating the distribution in this region, we can see that ∼45% of the zones are condensed. The 55% remaining is redistributed between endloops and some noncondensed internal zones. There are thus coexisting open and collapsed plectonemic regions on the same plasmid. This particular shape obtained at an effective diameter equal to the geometrical diameter of DNA has the very interesting property of presenting local zones with a high degree of tightness adjacent to open zones much more capable to interact with other biological molecules. This shape was already observed in cryoelectron microscopy (Bednar et al) and was also predicted in a statistical study of supercoiled DNA (Marko and Siggia, 1995).

Display large version of this figure
Figure 5
Two supercoiled DNA molecules of 2686 bp (pUC18) with a difference of linking number (ΔLk) equal to −13 in 20mM MgCl2 obtained with the introduction of an attractive potential into the simulation (A=0.86 kT). Note the coexistence of open and collapsed zones on the same molecule.
Display large version of this figure
Figure 6
Distribution of the superhelix diameter at 20mM MgCl2 of the pUC18 (ΔLk=−13) in aqueous solvent with the introduction of an attractive potential into the simulation (A=0.86 kT).

This cohabitation can be explained in the following way. It is known that the writhe is very sensitive to the effective diameter. With a low effective diameter the polymer can interwind more easily than with a larger one. For instance, we know already that a decrease of the effective diameter will give larger endloops (Schlick et al). Therefore, modeled molecules easily accommodate high writhe, which strongly reduces the effective gradient of torsional tension acting on supercoiled DNA molecules. Under such conditions there is no elastic energy gradient that could limit the entropy of noncollapsed regions. The coexistence of these two shapes is thus energetically favorable. It also confirms the fact that the reduction of conformational entropy is not completed with a degree of supercoiling equal to −0.05. A higher degree of supercoiling would be needed to condense the polymer entirely as also proposed by Marko and Siggia, 1995.

Interestingly, a long simulation time corresponding to 1ms relaxation of a molecule is needed until the initial circular configuration with ΔLk=−13 is equilibrated and reaches an interwound superhelix with interspersed collapsed and open regions. This time is ∼5 times longer than the equilibration time of molecules in simulations neglecting attractive interactions between DNA segments. Long equilibration time in simulations taking in account DNA-DNA attractive potential results from reduced possibilities of internal movements after some regions are laterally aligned. The “transversal entropy” is frozen with the combination of intermolecular attractive forces, and a high degree of supercoiling and reshaping of the molecule can thus only occur with “translational entropy,” e.g., reptation type of movements of segments along each other, and this is a slow process.

It was found by Hammermann et al by small-angle laser light scattering that the degree of the internal motion in supercoiled DNA begins to decrease from a high concentration of NaCl (1M). In accordance with our results, this phenomenon can be well-explained by the apparition of condensed zones facilitated by the high concentration of added monovalent salt. This idea is also reinforced by Monte Carlo simulations (Montoro and Abascal, 1995) and also by the application of the modified P-B model to DNA (Bhuiyan and Outhwaite, 1994) that showed an effect of “overneutralization” of the DNA charges in very concentrated monovalent electrolyte solutions. Experimental measurements of the internal motion of DNA molecules are a good way to study the existence of condensed zones in a supercoiled DNA, and further study by small-angle scattering, for instance at ∼20mM MgCl2, a concentration characterized by an effective diameter equal to the geometrical diameter of DNA, is needed. In a study by cryoelectron microscopy, Bednar et al had found that at elevated concentration of NaCl supercoiled DNA presented some zones of lateral contact. However, the exact value of salt concentration needed to induce the lateral contact was not known, because the concentration of salt is not well under control during the last second of the specimen preparation.


Simulated nonsupercoiled DNA molecules in 20mM MgCl2 in aqueous solvent do not condense despite taking into account attractive interactions

To check whether DNA supercoiling is required for the demonstrated above DNA collapse, we simulated behavior of nonsupercoiled DNA molecules in aqueous solutions containing 20mM MgCl2. The parameters of the simulation are kept the same as before except for the supercoiling density, which is now set to zero. Fig. 7 presents the evolution of the radius of gyration of a 1.5-kb circular nonsupercoiled molecule during 6ms in 20mM MgCl2 and water solvent. Note that we have chosen a smaller DNA molecule (1.5kb) in order to decrease computer time to have a sufficiently long run, as we know that site juxtaposition in DNA is a slow process (Jian et al). No supercoiling (ΔLk=0) was introduced, so that we can expect no condensation, as it is well-known from several experimental results. Indeed, simulated molecules do not seem to collapse in this studied time range because the radius of gyration does not present any sudden decrease. However, the distribution of the radius of gyration averaged over time presents a surplus of low values compared to a Gaussian distribution (see Fig. 8). A configuration obtained after 3.4ms and showing a low value of the radius of gyration is shown in Fig. 7. One can easily see a condensed zone on it. This intermediate state (the condensed zone has disappeared after a while) is a typical feature of a metastable state where entropic contributions act against pure enthalpic contributions. Thus the nonsupercoiled DNA passes through intermediate states where condensed zones are present, but released after a while by entropic effects. Lack of global condensation of nonsupercoiled DNA in 20mM MgCl2 is consistent with known literature data (Ma and Bloomfield, 1994).

Display large version of this figure
Figure 7
Evolution of the radius of gyration of a 1500 bp nonsupercoiled molecule during 6ms in 20mM MgCl2 starting from a circular conformation in aqueous environment without excess of torsional energy but accounting for the attractive potential (A=0.86 kT). A configuration with a locally collapsed zone was detected and shown after 3.4ms by observing significant decrease of the radius of gyration. No global aggregation was seen in this lapse of time.
Display large version of this figure
Figure 8
Distribution of the radius of gyration of a 1500 bp nonsupercoiled molecule in 20mM MgCl2 in aqueous solvent during 6ms of simulation accounting for the attractive potential (A=0.86 kT). The curve deviates from a gaussian distribution and shows regions with predominant small values of the radius of gyration.

Simulated nonsupercoiled DNA molecules in 20mM MgCl2 do condense in 50% methanol

It is known from an experimental study by light scattering (Post and Zimm, 1982) that DNA condensation occurs in two stages: first, an intramolecular collapse that can occur at already low concentrations of DNA, and second, an aggregation of many molecules at higher concentrations. The first stage is a fast reaction, in the millisecond time range (Porschke, 1984) and can be studied in this BD simulation.

The modification of the solvent with the addition of methanol is parametrized as follows. The introduction of 50% methanol decreases the dielectric constant by a factor of 0.7. This changes the two characteristic lengths of the electrostatic interactions, the Bjerrum length and the Debye length, which now take a value, respectively, of 1 and 1.05nm. We can note that the P-B formalism introduced in the theory that is based on a definition of an effective linear charge of DNA does not interfere with this change of the dielectric constant. In other words, the numerical factor multiplying the repulsive electrostatic potential (Eq. (9)) is only divided by a factor of 0.7 due to its inverse dependence on the dielectric constant.

As the attractive energy is predominantly a second-order electrostatic effect, it is in a first approximation proportional to the square of the Bjerrum length as also modeled by several authors (Oosawa, 1968,Barrat and Joanny, 1995,Gronbech-Jensen et al). Thus the effective attractive parameter A is taken as equal to 1.75 kbT. This parametrization gives an effective diameter below zero and we can thus expect to see condensation in this case. The other parameters are kept the same as for the previous simulations.

Fig. 9 presents the evolution of the radius of gyration of a simulated DNA molecule of 1500 bp in 50% methanol at 20mM MgCl2. We can see the sharp transition of the radius of gyration after ∼1.8ms. The process of condensation is relatively fast, with a decrease of the radius of gyration of ∼20 nm/0.2ms. Cooperative effects are clearly expected during this process, as we can see on the three configurations inserted in Fig. 9. The first configuration is obtained after the first juxtaposition of two subunits at 1.85ms. The second shows that during this simulation the condensation propagated along the DNA. The third shows further progress of DNA condensation. However, we have to note that upon decreasing the value of the radius of gyration below 30nm, some numerical instabilities appeared, preventing us from studies of the shape of the condensed molecules. These numerical problems are directly related to physical ones. For instance, a high number of helical forms were found during the condensation simulation. The helices of condensed DNA had a pitch of ∼2.5nm and a diameter of 10nm, which is only twice the length of a segment in our simulation. The study of these condensed zones, which are very similar to DNA around histones, requires a discretization better than 20 segments per statistical length, which would be very costly in computer time and also, below this value, molecular effects are clearly expected and put boundaries to our statistical approach. If we neglect these physical-numerical deficiencies, we can note that a compaction of a 500-nm polymer in a sphere of 13-nm radius was obtained, i.e., a very high degree of packaging.

Display large version of this figure
Figure 9
Evolution of the radius of gyration during 2.3ms of a 1500 bp nonsupercoiled molecule in 20mM MgCl2 in 50% methanol solvent. Precise parameters of the simulation can be found in the text. A sharp decrease of the radius of gyration is seen after the first contact between two subunits. The first inset at 1.85ms presents the conformation obtained when this first contact between two monomers occurs. Aggregation propagates from this point along the molecule as shown by the two pictures at time 2.1ms and 2.25ms, respectively.

Simulation of a nonsupercoiled DNA molecule in 20mM MgCl2 confined in a box

It is known that a very high concentration of DNA can lead to an aggregation of the molecule even in a very small amount of NaCl (Wissenburg et al). Furthermore, we have shown that a reduction of entropy caused by a high degree of supercoiling in conditions where the second virial coefficient gives an effective diameter equal to the geometrical diameter can lead to a local intermolecular collapse. Therefore, we study the effect of confining a DNA molecule in a cube whose size is variable. For this, a 1500-bp circular unsupercoiled molecule is plunged initially into a cube with side size equal to the diameter of the fully unfolded circular DNA molecule (see first inset in Fig. 10). Then we linearly and isotropically decrease the size of the box. The potential of confinement is introduced to reproduce an infinite potential outside the cube. To avoid numerical divergences it is evaluated as an exponential with a decay rate much smaller than each dimension in the simulation (presently 0.1nm). Furthermore, in order to allow a possible internal condensation, the rate of decreasing the cube size has to be much slower than the speed of reshaping of the molecule. However, for practical reasons of required computer time we chose a compromise with a decreasing rate of 20nm per ms of simulation. Fig. 10 presents both the evolutions of the radius of gyration and of the percentage of collapsed zones. We can see that at the beginning of the simulation the fluctuations of the radius of gyration are not limited by confining effects. Behavior of the DNA molecule is comparable to a random coil, and no collapsed zones are found in the molecule. From ∼4.3ms another stage appears in the simulation (see second inset in Fig. 10). Confinement restricts possible movements of the molecule and some local collapsed zones are clearly visible. In further decreasing the size of the cube side, the molecule is pushed against the walls until it reaches an aggregated conformation with a striking resemblance with a torus (see third inset at 6.9ms in Fig. 10). In this simulation, two different mechanisms of condensation have to be taken in account: first an internal local collapse, which can be obtained by the combination of small attractive forces and a reduction of conformational entropy; and second, a global aggregation of the molecule due to the direct action of confining the molecule. As both these effects reduce the distance between monomers, it is difficult in this simulation to disentangle them. However, we have to note that no clear transition of the radius of gyration was seen even when some zones of lateral contact are clearly visible from ∼4.3ms, showing that the final aggregated form is only obtained with the help of external pressure of confinement.

Display large version of this figure
Figure 10
Evolutions of the radius of gyration and of the percentage of collapsed zone of a 1500 bp nonsupercoiled (ΔLk=0) simulated molecule in 20mM MgCl2 in aqueous solvent in a shrinking cube accounting for the attractive forces (A=0.86 kT). Half of the cubes side is also shown for comparison. Three conformations are shown at respectively 0, 4.3, and 6.9ms. Proportions are respected in these pictures with a diameter of the molecule corresponding to 2nm. We define a collapsed zone as a zone where the distance between opposing segments is smaller than 4nm.

During the simulation of the supercoiled molecule in 20mM MgCl2 we have shown that ∼45% of zones along one molecule were collapsed with a degree of supercoiling equal to −0.05. This fraction of collapsed zones can be compared with the present simulation. Such a percentage is obtained after ∼6.3ms of simulation with a corresponding concentration of DNA of ∼160 g/l. Such a high value confirms that supercoiling is a powerful tool for compacting DNA.


General comparison with recent experimental results

We have demonstrated here using BD simulations that naturally supercoiled DNA molecules adopt partially collapsed structure in aqueous solutions where concentrations of divalent Mg2+ cations approach and exceed 20mM. Both the supercoiling and the relatively weak attractive DNA-DNA interactions that are known to operate under these conditions are needed for the collapse. It is an accepted fact that counterions with a valency higher than two are needed for a condensation of linear DNA in aqueous solutions, and our simulations of long circular unsupercoiled DNA in 20mM Mg2+ correctly demonstrated that weak attractive DNA-DNA interactions acting under these conditions are unable to compensate for the entropy loss needed for the condensation of linear DNA molecules. Of course, simulations are only simulations, and we are aware that the outcome of the simulation can be affected by the uncertainty of some entered parameters and by substantial simplifications of interactions and dynamics of the system that were needed to facilitate computations. Therefore, a physical experimental verification of our results is required. Electron microscopy provided many images of collapsed DNA molecules under conditions where the imaged linear DNA molecules do not collapse, but it was frequently argued that electron microscopy is not a good method for concluding the state of DNA in solution. Even the results obtained by cryoelectron microscopy, where one can observe rapidly vitrified specimens with supercoiled DNA molecules immersed in a buffer of interest, are called in question because the salt concentration is not under control. Therefore, solution methods are needed to study DNA-DNA attraction in the presence of divalent cations. Shaw and Wang, in their studies of DNA knotting upon end-closure of long linear DNA molecules, were able to precisely calibrate the effective diameter of DNA in buffers with increasing concentrations of Mg2+. The presence of DNA-DNA attraction is manifested by bringing down the effective diameter of DNA below the 2nm value, i.e., the value of the geometric diameter of DNA. Also, studies of supercoiled DNA performed by Rybenkov et al demonstrated that under these conditions the effective diameter of DNA becomes smaller than the geometric diameter of DNA. While the value of the effective diameter smaller than this of the geometric diameter indicates the presence of attractive DNA-DNA interactions, it is still an open question whether these attractive interactions are strong enough to induce an internal collapse of supercoiled plasmids. It turns out that it is a difficult task to conclude from standard physicochemical experimentation whether there is an internal DNA-DNA association within supercoiled DNA molecules.

Two recent papers by Rybenkov et al,Rybenkov et al addressed the question of the collapse of supercoiled DNA by comparing experimental results of indirect experimental methods with the simulated behavior of supercoiled DNA molecules for which no intersegmental attraction was entered into simulation. These authors observed a good agreement between the experimentally measured sedimentation coefficient of the supercoiled DNA molecules in high salt solution and the calculated sedimentation coefficient for the simulated DNA molecules. Since the models used for the simulation of DNA behavior operated only with a simple hard-core repulsive potential, the obtained simulated configurations were free from the regions of collapse. The fact that the calculated sedimentation coefficient for the simulated molecules without the region of collapse matched the sedimentation coefficient of the supercoiled DNA molecules in high salt solution was used by Rybenkov et al as an argument that there is no collapse of the supercoiled DNA in high salt solutions. However, we do not know if the fit between experimental and simulated values would have been worse if the attractive potential leading to the collapse of the supercoiled DNA was introduced into the simulations. It is likely that the sedimentation coefficient is a rather robust property that can be practically independent from the actual distance between the opposing segments of the interwound superhelices. In fact, Rybenkov et al clearly indicated that this is the case by writing that the sedimentation coefficient is “much less sensitive, if at all, to the superhelix diameter.” Thus it is very surprising that Rybenkov et al reject the possibility of DNA collapse when a local collapse may be fully compatible with results obtained by them. We emphasize the fact that we have to distinguish between a local collapse, which can be difficult to observe by indirect methods, and a global aggregation easily visible experimentally because of the domination of pure attractive enthalpic effects over entropic ones. In their second paper, Rybenkov et al measured the efficiency of catenanes formation upon circularization of linear duplex DNA in the presence of supercoiled DNA. They have observed that the frequency of threading a linear duplex between the opposing supercoiled DNA increased with increasing concentration of salts. These results were interpreted by the authors as the proof that there is no collapse, as otherwise it would be difficult to spontaneously thread linear duplex between the opposing segments of supercoiled DNA that are supposedly in continuous direct contact with each other. However, Rybenkov et al did not consider that lateral association of the opposing segments does not imply that these segments are strongly bound to each other. In fact, micrographs obtained by the cryo-EM method (Bednar et al), a statistical study of supercoiled DNA (Marko and Siggia, 1995), and our present simulations demonstrate that the regions of collapse are interspersed with open regions having a relatively big cross-sectional area. In addition, terminal loops in collapsed supercoiled DNA are bigger than in noncollapsed DNA because it is also for small effective diameter (Schlick et al) and thus there is plenty of space for spontaneous threading of linear DNA through “collapsed” supercoiled DNA leading to formation of catenanes. This space is increased by the fact that at these conditions the effective diameter of DNA is very small and the loops are many times larger that the diameter of the threading molecule. Furthermore, catenation efficiency, studied by Rybenkov et al, is also insensitive to conditions where it is known that an intermolecular aggregation occurs (3.5mM spermidine). Taking into account the above-mentioned results we can say that increased probability of catenation between supercoiled and linear DNA under the condition of substantial charge neutralization is fully consistent with the collapse observed in our simulations.

A recent study by light scattering (Hammermann et al) allowed precise measurement of the internal dynamics of superhelical DNAs. They showed that from ∼1M NaCl the internal motion of the DNA molecule presented a sudden decrease, and interpreted this result as a sign that a compaction of the overall DNA structure occurs at a high concentration of monovalent salts. Since divalent counterions have ∼50-fold higher screening potential than monovalent ones (Shaw and Wang, 1993), this interpretation is consistent with the results presented here. A similar study in divalent solutions like MgCl2, which would also present a sudden decrease of the internal motion from ∼20mM, could reinforce the idea of local collapse as developed in this paper.



Conclusions

By using the BD simulation approach incorporating a simple model of DNA-DNA attractive forces, we have shown that it is possible that a high reduction of entropy caused either by supercoiling or by confining the molecule in a small box can facilitate a local intramolecular condensation of the DNA even in 20mM MgCl2, i.e., close to physiological concentrations. As at the θ-point, the second virial coefficient has the interesting property of being characterized at 20mM MgCl2 by a cancellation of the repulsive electrostatic forces by attractive ones. The BD simulation is a very useful tool in understanding all forms of condensation of polymers. The simulation takes into account geometrical and topological features of the polymer, repulsive and attractive electrostatic interactions between monomers, hydrodynamical interactions, and entropic effects, which are all responsible directly or indirectly for the polymer condensation. However, limits are given by the introduction of precise independent parameters and by the fact that at small dimensions molecular effects have to be taken into consideration. The introduction of the MPB equation should be the next step and will allow us to study the shape of condensed polymer systems and their possible phases or their loss of isotropy upon condensation. Because we know that the attractive potential is very sensitive to the steric size of the counterion, this formalism could also explain the still poorly understood fact that such divalent counterions as Mn2+ can lead to an intermolecular aggregation of supercoiled DNA molecules, although Mg2+ cannot (Ma and Bloomfield, 1994).


Acknowledgments

This work was supported by Swiss National Science Foundation Grant 31-42158.94 and by Human Frontier Science Program Grant RG 335/97.

References

Allison et al., 1989 Allison, S., Austin, R., and Hogan, M. (1989). Bending and twisting dynamics of short linear DNAs. Analysis of the triplet anisotropy decay of a 209 base pair fragment by Brownian simulation. J. Chem. Phys. 90, 3843–3854. CrossRef | PubMed

Barrat and Joanny, 1995 Barrat, J.-L., and Joanny, J.-F. (1995). Theory of polyelectrolytes solutions. In Prigogne, I., Rice, S., eds. Advances in Chemical Physics 94, (New York: John Wiley and Sons), pp. 1–66. PubMed

Bednar et al., 1994 Bednar, J., Furrer, P., Stasiak, A., Dubochet, J., Egelman, E.H., and Bates, A.D. (1994). The twist, writhe and overall shape of supercoiled DNA change during counterion-induced transition from a loosely to a tightly interwound superhelix. J. Mol. Biol. 235, 825–847. CrossRef | PubMed

Bhuiyan and Outhwaite, 1994 Bhuiyan, L.B., and Outhwaite, C.W. (1994). The cylindrical electric double layer in the modified Poisson-Boltzmann theory. Philosophical Magazine B 69, 1051–1058. PubMed

Bloomfield, 1996 Bloomfield, V.A. (1996). DNA Condensation. Curr. Opin. Struct. Biol. 6, 334–341. CrossRef | PubMed

Bloomfield, 1997 Bloomfield, V.A. (1997). DNA Condensation by multivalent cations. Biopolymers/Nucleic Acid Sci. 44, 269–282. PubMed

Bloomfield et al., 1980 Bloomfield, V.A., Wilson, R.W., and Rau, D.C. (1980). Polyelectrolyte effects in DNA condensation by polyamines. Biophys. Chem. 11, 339–343. CrossRef | PubMed

Brenner and McQuarrie, 1973 Brenner, S.L., and McQuarrie, D.A. (1973). Force balances in systems of cylindrical polyelectrolytes. Biophys. J. 13, 301–331. Abstract | | CrossRef | PubMed

Chirico and Langowski, 1994 Chirico, G., and Langowski, J. (1994). Kinetics of DNA supercoiling studied by Brownian dynamics simulation. Biopolymers 34, 415–433. CrossRef | PubMed

Clement et al., 1973 Clement, R.M., Sturm, J., and Daune, M.P. (1973). Interaction of metallic cations with DNA VI. Specific binding of Mg2+ and Mn2+. Biopolymers 12, 405–421. CrossRef | PubMed

Das et al., 1997 Das, T., Bratko, D., Bhuiyan, L.B., and Outhwaite, C.W. (1997). Polyelectrolyte solutions containing mixed valency ions in the cell model: a simulation and modified Poisson-Boltzmann study. J. Chem. Phys. 107, 9197–9207. CrossRef | PubMed

Duguid et al., 1993 Duguid, J., Bloomfield, V.A., Benevides, J., and Thomas, G.J. (1993). Raman spectroscopy of DNA-metal complexes. I. Interactions and conformational effects of the divalent cations: Mg, Ca, Sr, Ba, Mn, Co, Ni, Cu, Pd and Cd. Biophys. J. 65, 1916–1928. Abstract | | CrossRef | PubMed

Ermark and McCammon, 1978 Ermark, D.L., and McCammon, J.A. (1978). Brownian dynamics with hydrodynamic interactions. J. Chem. Phys. 69, 1352–1360. CrossRef | PubMed

Felgner and Ringold, 1989 Felgner, P.L., and Ringold, G.M. (1989). Cationic liposome mediated transfection. Nature 337, 387–388. CrossRef | PubMed

Gavryushov and Zielenkiewicz, 1998 Gavryushov, S., and Zielenkiewicz, P. (1998). Electrostatic potential of B-DNA: effect of interionic correlations. Biophys. J. 75, 2732–2742. Abstract | Full Text | PDF (191 kb) | CrossRef | PubMed

Gronbech-Jensen et al., 1998 Gronbech-Jensen, N., Beardmore, K.M., and Pincus, P. (1998). Interactions between charged spheres in divalent counterion solution. Physica A 261, 74–81.