Quantifying the Thermodynamic Consequences of Cation (M2+, M+)

Quantitative interpretations of the large Coulombic effects of changes in cation concentrations on processes involving polyanionic DNA require accurat...
0 downloads 0 Views 183KB Size
J. Phys. Chem. B 1999, 103, 3489-3504

3489

Quantifying the Thermodynamic Consequences of Cation (M2+, M+) Accumulation and Anion (X-) Exclusion in Mixed Salt Solutions of Polyanionic DNA Using Monte Carlo and Poisson-Boltzmann Calculations of Ion-Polyion Preferential Interaction Coefficients Haihong Ni,† Charles F. Anderson,† and M. Thomas Record, Jr.*,†,‡ Departments of Chemistry and Biochemistry, UniVersity of WisconsinsMadison, Madison, Wisconsin 53706 ReceiVed: NoVember 10, 1998; In Final Form: January 26, 1999

Quantitative interpretations of the large Coulombic effects of changes in cation concentrations on processes involving polyanionic DNA require accurate theoretical descriptions of the thermodynamic consequences of cation-DNA interactions. In the present study, the thermodynamic consequences of accumulation of divalent and/or univalent cations in the vicinity of polyionic double-stranded DNA and of exclusion of univalent co-ions are characterized by ion-polyion preferential interaction coefficients Γi (i ) 2+, +, or -). These are calculated using integrals of ion distributions generated from either canonical Monte Carlo (CMC) simulations or numerical solutions of the cylindrical Poisson-Boltzmann (PB) equation, for the same minimally parameterized cylindrical cell model over experimentally relevant ranges and ratios of the uni- and/or divalent cations. For solutions containing both types of cations, trends in ΓiMC and in ΓiPB with changes in the absolute and relative values of the divalent and univalent cation concentrations are examined and compared with trends calculated for solutions containing only one type of cation. Differences between ΓiMC and ΓiPB are quantified and related to differences between MC and PB predictions of the extents of local cation accumulation within 3 Å of the polyion surface. Discrepancies between Γ2+MC and Γ2+PB, and between Γ-MC and Γ-PB, are significant whether or not univalent cations are present, but the difference ∆Γ2+ ) Γ2+MC - Γ2+PB is relatively insensitive to changes in the concentration of salt (excess 1:1). Therefore, PB calculations may provide a satisfactory alternative to more computationally demanding MC simulations as a basis for analyzing the saltconcentration dependences of Γ2+ and of the closely related measurable thermodynamic properties that reflect the importance of Coulombic interactions.

I. Introduction The present study is motivated by the need for theoretically rigorous calculations of the thermodynamic consequences of Coulombic interactions of salt ions, such as Na+, Mg2+, and Cl- with polyanionic DNA in mixed salt solutions, for use in the analysis of effects of changes in concentration of these ions on the equilibria and kinetics of nucleic acid processes. In Escherichia coli and other cells, ligand-nucleic acid interactions and the formation of ordered nucleic acid structures occur in mixed salt polyelectrolyte solutions containing a relatively high concentration (g10-1 M) of unneutralized nucleic acid polyanion charge and moderate concentrations (10-3-10-2 M) of free Mg2+ and oligocationic “polyamines” (especially putrescine (2+) in E. coli) as well as an osmotically variable concentration of K+ (>10-1 M).1 Physiological concentrations (mM) of Mg2+ or polyamines are required to form the ribosome and other folded RNA structures and, in many cases, to obtain functional protein-DNA complexes (e.g., complexes with restriction enzymes, other nucleases, and polymerases). For this reason, and/or because of the requirements of the binding assay, protein-nucleic acid interactions and nucleic acid conformational changes are often investigated in mixed salt solutions * To whom correspondence should be addressed: Department of Biochemistry, 433 Babcock Drive, University of WisconsinsMadison, Madison, WI 53706. Phone: (608) 262-5332. Fax: (608) 262-3453. E-mail: [email protected]. † Department of Chemistry. ‡ Department of Biochemistry.

(containing both 1:1 and 2:1 electrolytes, such as NaCl and MgCl2, respectively). Large nonspecific (Coulombic) effects of Mg2+ are generally observed on nucleic acid conformational transitions and proteinnucleic acid binding interactions, which also may manifest specific effects of Mg2+ (as indicated above). An increase in the concentration of Mg2+ stabilizes helical states of nucleic acids and destabilizes both site-specific and nonspecific proteinDNA complexes. In the absence of 1:1 salt, these Coulombic effects of Mg2+ require a polyelectrolyte analysis, but in excess 1:1 salt, they have been analyzed successfully using a mass action equilibrium model with a Mg2+-DNA binding constant that exhibits a power-law dependence on the concentration of excess NaCl.2-4 A more fundamental thermodynamic analysis of these large Coulombic effects of cations on processes involving nucleic acid polyanions requires an accurate theoretical description of the thermodynamic consequences of changes in the local distributions of salt ions in response to changes in bulk salt concentrations. Molecular and thermodynamic consequences of ion-polyion interactions can be quantified, analyzed, and interpreted in terms of preferential interaction coefficients, which also provide a novel and general method of characterizing the contributions of polyelectrolyte, Hofmeister, and osmotic effects to the observed salt-concentration dependences of all processes involving biopolymers (proteins and/or nucleic acids).1 Preferential interaction coefficients have been defined in various (closely related) forms that are appropriate for different kinds of systems,

10.1021/jp984380a CCC: $18.00 © 1999 American Chemical Society Published on Web 04/09/1999

3490 J. Phys. Chem. B, Vol. 103, No. 17, 1999 experimental constraints, and methods of measurement or calculation. Alternative representations of preferential interaction coefficients, relationships among them, and connections with quantities directly accessible to experimental measurement have been discussed in detail,5-8 primarily with regard to the interactions of uncharged solutes or electroneutral components comprised of dissociated ionic species. When defined as derivatives of the molal concentrations of electrolyte and polyelectrolyte components, preferential interaction coefficients form a basis for the thermodynamic analysis of observable effects of excess salt concentration on processes at constant temperature and pressure, such as DNA-oligocation binding.9 (The appendix of ref 9, for example, relates this type of a preferential interaction coefficient to the corresponding Donnan coefficient determined by membrane equilibrium.) To characterize counterion accumulation and co-ion exclusion arising from salt ion-polyion interactions in solutions containing one type of electrolyte and one type of polyelectrolyte component, preferential interaction coefficients for individual salt ions were defined10 in the context of the Donnan ion distribution across a membrane impermeable to polyions. Equivalent thermodynamic functions, called ion-polyion preferential interaction coefficients and designated Γi (i ) +, 2+, or -) in the present paper, are evaluated herein by integrating the radial distributions calculated for univalent and/or divalent cations and for univalent co-ions surrounding polyanionic DNA in solutions that contain one or two added salts (e.g., NaCl and/or MgCl2). These theoretical predictions of the salt ion distribution functions are generated both from canonical Monte Carlo (CMC) simulations and from numerical solutions of the Poisson-Boltzmann (PB) equation for the standard cylindrical cell model of polyanionic DNA in aqueous solution with added salt(s). For polyelectrolyte solutions that contain both 1:1 and 2:1 salts, the cylindrical PB equation was used previously11,12 to calculate a function that is related closely to Γ2+ for the purpose of analyzing experimentally determined quantities. Previously reported MC calculations of preferential interaction coefficients, defined in terms of electroneutral components rather than individual ionic species, have been based on the grand canonical ensemble and pertain to polyelectrolyte solutions containing only a single type of a 1:1 electrolyte.13-16 Our CMC results provide, for the first time, tests of the accuracy of PB-based predictions of the ion-polyion preferential interaction coefficients that quantify the thermodynamic consequences of the interactions of divalent cations with a DNA polyanion in a 2:1 salt solution, in either the presence or absence of 1:1 salt. The following section (II) summarizes some essential aspects that enter into the calculation of ion-polyion preferential interaction coefficients by either of the theoretical methods employed herein. Section III specifies the model assumptions common to both of these methods and explains how input from CMC simulations or numerical solutions of the cylindrical PB equation was obtained and used to calculate ion-polyion interaction coefficients. Results over a representative range of salt concentrations are presented and discussed in section IV, and their physical implications are placed in the broader context of previous work in section V. Applications of these ionpolyion preferential interaction coefficients in thermodynamic analyses of cation binding to DNA and of effects of salt concentration on the equilibria and kinetics of nucleic acid processes are made in subsequent papers (Ni et al., in preparation).

Ni et al. II. Calculating Ion-Polyion Preferential Interaction Coefficients as Integrals over Ion Distributions Taking into account the relevant characteristics of the systems and conditions investigated herein, we consider the physical basis for calculating an ion-polyion preferential interaction coefficient, Γi, as the difference between an integral over the nonuniform salt ion distribution around an isolated cylindrical polyion and an integral over the uniform ion distribution in a corresponding reference solution, where interactions involving the polyion are absent. To implement our calculations of Γi, the requisite ion distributions were obtained either from CMC simulations (as described in section IIIB) or from numerical solutions of the cylindrical PB equation (as described in section IIIC). When these alternative sources of input are used, different factors affect the calculation of ΓiPB or ΓiMC, as explained in section IIID. In all of the solutions for which we have calculated Γi, the concentration of the polyionic species is sufficiently dilute so that (in the statistically significant configurations) virtually no salt ions interact with more than one polyion and most salt ions interact, in effect, only with other salt ions. Consequently, the ion distributions around each polyion are unaffected by interactions with any other polyion. In our MC and PB calculations, the polyion is modeled as a uniformly and continuously charged cylinder that is long enough so that the spatial distributions of the salt ions surrounding the polyion can be calculated, per polyion charge, without considering Coulombic end effects.15-17 Because the terminal regions of the polyion where these effects are manifested constitute a negligible portion of its total length, each ion distribution can be represented as a function of one variable: the radial coordinate r, zeroed at the polyion axis. Specifically, for each type of ion i the local density, Ci(r), in units of molarity, is represented as a function of r, in units of angstroms, by defining 2π(10-27NA)Ci(r) r dr (where NA is Avogadro’s number) as the differential increment of ions whose centers fall within the annular shell extending from r to r + dr. The polyion and salt ions are modeled as impenetrable (“hard”) particles. Therefore, Ci(r) ) 0 for 0 e r < a, where a is the distance of closest approach of the center of a small ion to the polyion axis. In accordance with the “restricted” primitive model, each type of salt ion has the same radius. (The radii assigned to the polyion and to the small ions, together with the values of all of the other parameters needed to specify the model, are given in section IIIA.) For r > a, Coulombic attractions and repulsions, mediated by the dielectric properties of pure water, produce gradients in the local densities of both positively and negatively charged salt ions surrounding an isolated polyion. For the systems and conditions investigated herein, the steepness of each of these gradients is determined primarily by the average (projected) axial density of structural charges on the rodlike polyion and by the concentration(s) of salt(s). The concentration of the polyionic species is sufficiently dilute so that far from any polyion the systematic trends in the Ci(r) produced by screened Coulombic interactions with the polyion charges become negligible compared with fluctuations in the local densities of the salt ions, and each of the Ci(r) has reached the corresponding spatially uniform value, Cibulk. In a reference system whose thermodynamic state is determined entirely by the properties of the bulk regions of the polyelectrolyte-salt solution, Ci(r) ) Cibulk at every r g 0. Consequently, in systems where ion-polyion interactions are operative, the magnitude of the difference Ci(r) - Cibulk is a spatially dependent measure of effects that arise from Coulombic attractions and repulsions and from the volume excluded to salt

Mixed-Salt Solutions of Polyanionic DNA

J. Phys. Chem. B, Vol. 103, No. 17, 1999 3491

ions by the polyion. For each type of ion i, the total extent of accumulation or exclusion produced by these interactions can be expressed as the ion-polyion preferential interaction coefficient, Γi, obtained by integrating Ci(r) - Cibulk over (at least) the entire annular volume of solution where these concentrations are significantly different. The radius enclosing such an annular volume herein is designated Rd. The following mathematically equivalent alternative integral representations are useful from different theoretical or methodological standpoints:

∫0R r(Ci(r) - Cibulk) dr

Γi ) 2π(10-27NA)b

d

(1a)

∫aR r(Ci(r) - Cibulk)dr -

) 2π(10-27NA)b

d

π(10-27NA)bCibulk a2 (1b)

∫0R rCi(r) dr -

) 2π(10-27NA)b

d

π(10-27NA)bCibulk Rd2 (1c) ) {Ni}d - {Nibulk}d

(1d)

In these equations, b (in units of angstroms) is the average axial spacing between monomeric charges on the cylindrical polyion. Thus, Γi is expressed in nondimensional form per polyion monomer charge. Equation 1a indicates that Γi can be calculated by direct numerical integration of ion distributions obtained from numerical solutions of the cylindrical PB equation (as explained in section IIIC) or from any other method capable of predicting the spatial dependence of Ci(r). The lower limit on this integration is set at r ) 0, to include the contribution that arises solely from the volume excluded to salt ions by the polyion monomer. This contribution is shown explicitly as the second term on the right of eq 1b. (The first term in this equation also depends on the structural radius of the polyion, except in solutions so dilute that the Γi have practically attained their limiting values, which depend on b, but not a.)18 With the discretized ion distributions obtained from MC simulations (as described in section IIIB), values of ΓiMC can be calculated on the basis of eq 1c. This equation is expressed more compactly as eq 1d by introducing the following symbols: {Ni}d, for the averaged number of ions i per polyion monomer that are located in the annular volume per polyion monomer extending in the radial direction to Rd; {Nibulk}d, for the number of ions i that would be located in the volume πbRd2 if all types of interactions involving the polyion were absent (i.e., if it were completely penetrable, with no charge). As calculated from any of the expressions in eq 1, the magnitude of Γi is unaffected by an increase in the upper limit of the integral as long as Ci(Rd) remains equal to Cibulk. This equality holds for any value of Rd that exceeds the effective range of the potential of mean force characterizing the interactions of an isolated polyion with any of the salt ions, provided also that Rd is less than half the average distance between nearest-neighbor polyions. (Otherwise, the neglect of interpolyion interactions in the model used to calculate the Ci(r) may not be justified.) For an infinitely long polyion cylinder, at infinite dilution, Rd could be extended to infinity, and eq 1a could be used to evaluate Γi as a kind of second virial coefficient19,20 by introducing numerical solutions of the cylindrical PB equation and an approximate analytic expression for the difference Ci(r) - Cibulk, valid at large values of r where the PB equation can be linearized. Instead of the infinite cell (IC) boundary conditions employed in Stigter’s classic study,20

our calculations of ΓiPB are based on finite cell (FC) boundary conditions that match exactly those used in our corresponding calculations of ΓiMC. In section IIID, we indicate how values of Rd were chosen to implement our use of the appropriate form of eq 1 to calculate ΓiMC and ΓiPB for model systems enclosed within finite electroneutral cells. The calculations reported herein pertain only to solutions that (in the statistically significant configurations) contain “bulk” regions far from any polyion that are electroneutral to a sufficient approximation (∑iziCibulk ) 0). Because the solution as a whole also is electroneutral, each polyion is surrounded by a (virtually) electroneutral annular volume that encloses the entire range over which ion-polyion interactions have a significant effect on any of the Ci(r). Within any finite electroneutral polyion-centered annular volume, where some or all of the Ci(r) are nonuniform, the average concentration of each type of salt ion differs from its corresponding bulk value. So, in accordance with eq 1c, Γi can be expressed as the quotient (C h i - Cibulk)/Cu. Here C h i is the average concentration of ions i located in the volume of solution that encloses exactly one polyion monomer, and the molar concentration of these monomers is Cu. If Cu is sufficiently dilute, a change in Cu while each of the Cibulk remains (virtually) constant has no effect on any of the Γi. Such a change in Cu would cause C h i to vary as a linear function of Cu with a slope Γi whose magnitude, at constant T and P, is determined only by the concentration(s) of the salt component(s) in the electroneutral bulk regions. On the basis of the foregoing considerations, Γi, calculated as an integral according to any of eq 1, is an ion-polyion preferential interaction coefficient whose magnitude reflects the extent of counterion accumulation or coion exclusion. In polyelectrolyte solutions containing one or more salt components, the set of Γi calculated for the individual salt ions are interrelated because both the solution as a whole and (to a sufficient approximation) the bulk regions far from any polyion are electroneutral. If the solution contains only a 1:1 salt and a polyion with univalent counterions:

Γ- ) Γ1:1 ) Γ+ - 1

(2)

Here, Γ1:1 denotes a preferential interaction coefficient defined in terms of the molar concentrations of the 1:1 salt component and the polyelectrolyte (monomer) component. With Γ2:1 defined analogously for a solution containing only 2:1 salt and a polyanion with divalent counterions, the relationship corresponding to eq 2 is:

Γ- ) 2Γ2:1 ) 2Γ2+ - 1

(3)

For a polyelectrolyte solution that contains univalent and divalent cations, both at concentrations exceeding Cu, counterions of either type can be assigned to the electroneutral polyelectrolyte component. Regardless of how this assignment is made, preferential interaction coefficients pertaining to the individual salt ions are interrelated by the condition of electroneutrality:

Γ- ) 2Γ2+ + Γ+ - 1

(4)

In solutions containing more than one type of salt component, relationships of the Γi to measurable quantities are determined by constraints characteristic of the experimental situation. The connection between experimentally accessible quantities and ion-polyion preferential interaction coefficients in DNA solutions containing one or two salt components will be addressed in a subsequent paper (Ni et al., in preparation).

3492 J. Phys. Chem. B, Vol. 103, No. 17, 1999

Ni et al.

TABLE 1: CMC Simulation Parameters for Solutions Containing a Model DNA Polyion (b ) 1.7 Å) and 2:1 or 1:1 Salt number of ionic species MC run number

salt

C h salt (mM)

cations

anions

polyion charges

approach distance a (Å)

cell radius Rc (Å)

cutoff radius Rd (Å)

fraction of accepted moves

1 1A 2 3 3A 4 5 6 6A 7 8 8A 9 10

2:1 2:1 2:1 2:1 2:1 2:1 2:1 1:1 1:1 1:1 1:1 1:1 1:1 1:1

13.4 13.4 72.0 113 113 169 239 16.5 16.5 72.0 113 113 169 239

130 130 180 200 200 300 210 160 160 180 200 200 300 210

200 200 300 340 340 540 390 100 100 120 140 140 240 180

60 60 60 60 60 60 30 60 60 60 60 60 60 30

10 13 10 10 13 10 10 10 13 10 10 13 10 10

229 229 118 100 100 100 100 229 229 118 100 100 100 100

75 75 60 55 55 45 40 160 160 70 65 65 60 50

0.526 0.553 0.577 0.535 0.527 0.558 0.497 0.591 0.588 0.565 0.583 0.571 0.603 0.571

Particularly for the ΓiMC, eqs 2-4 provide useful checks on the accuracy with which these coefficients have been evaluated by integrating the individual ion distributions according to eq 1c. The insensitivity of individual coefficients to the upper bound assigned to the integral (Rd) also can be checked by showing that the set of Γi (as calculated by either theoretical method) satisfy the appropriate condition of electroneutrality. III. Model and Methods (A) Specification of the Model. The univalent or divalent cations and univalent anions, representing hydrated monatomic ions Na+, Mg2+, and Cl- (for example), are modeled as hard spheres with equal radii: σ ) 2 or 3 Å. Double-stranded DNA is modeled as a cylinder with a continuous and uniform axial charge spacing b ) 1.7 Å and radius ao ) 8 or 10 Å. Corresponding to the larger and the smaller choices of σ and ao, MC and PB calculations were performed for two values of a ) ao + σ: 10 and 13 Å. The radii assigned to the polyion and small ions may be considered to include some water of hydration, but solvent water is modeled as a structureless continuum. The dielectric properties of the system are modeled as spatially uniform, locally and macroscopically concentrationinvariant, and equivalent to those of pure water (with dielectric constant  ) 78.4 at 25 °C). Calculations using a twodimensional form of the PB equation indicate that reducing the dielectric constant of the polyion interior (to 2, while retaining  = 80 for the surrounding ionic solution) produces relatively minor effects on ion distributions,17 which should become negligible when these are integrated to obtain ion-polyion preferential interaction coefficients. (B) Calculation of Ion Distributions from Canonical Monte Carlo Simulations. For each salt-polyelectrolyte solution modeled with the structural parameters specified above, the MC simulations from which we calculated ion-polyion radial distributions were implemented using canonical ensemble constraints in substantially the same way as we have described in detail previously.21 More general background on the MC method can be found in the appropriate section of the book by Allen and Tildesley.22 CMC simulations yield information about the distribution of each type of salt ion around an isolated polyion in the form of the average number of ions, Ni(r), located within each of a discrete set of annular volumes whose distance from the polyion axis is indicated by r. Specifically, in our simulations Ni(r) refers to the average number of ions whose centers lie within an annular volume with radial boundaries extending from (r - 0.2) to (r + 0.2) Å. The radial midpoint of the volume immediately adjacent to the polyion surface is

designated s ≡ a + 0.2. The concentration of ions i located in this volume is called a surface concentration and symbolized by Ci(s). At the temperature T ) 25 °C, the cylindrical MC cell contains a fixed total number Ni of each type of mobile salt ion (i ) +, 2+, or -) and a polyion segment, comprised of Np equally spaced discrete univalent charges, that is fixed along the axis of the cell. For each type of ion, the average concentration C h i is defined as Ni/Vc, where the cylindrical cell volume Vc includes the space occupied by the polyion segment. For each of the simulations in the present study, Tables 1 and 2 specify the values of the concentrations of the electroneutral h 2:1, all in the range 0.01-0.25 M, salt components C h 1:1 and C together with the number of salt ions of each type Ni and the dimensions (in angstroms) that determine the MC cell volume Vc ) 2πhcRc2. The axial length of the cell is 2hc ) 1.7Np, and Rc is the maximum radial distance from the polyion axis that can be occupied by the center of a small ion. According to the standard cylindrical cell expression,23 Rc is related to the polyion monomer molar concentration, Cu ) 1027Np/NAVc, by

Rc ) x1027/NAπbCu

(5)

In section IIID, we consider further the physical significance of Rc as an upper bound on the value of Rd in the integral used to evaluate ΓiPB or ΓiMC with the appropriate form of eq 1. A practical upper bound on the number of mobile salt ions enclosed within Vc is dictated by the amount of computational time required to maintain an accounting of the continually changing contributions to the configurational energy made by all energetically significant pairwise interactions experienced by each ion that is moved at each step during the simulation. In all simulations reported herein, both the axial and radial dimensions of Vc were large enough so that no further increases in Rc or in hc had any detectable effect on any of the Ni(r). The finding that ΓiMC does not change with decreasing Cu (outside statistical uncertainty) implies that this coefficient has the value characteristic of the limit Cu f 0. In the axial direction, longrange interactions between each ion within the cell and charges (either fixed on the polyion or on mobile salt ions) outside the cell are taken into account by supplementing the minimum image criterion22 with a self-consistently calculated external potential as described by Mills et al.21 In all of our simulations, hc was large enough to ensure that within each of the discrete annular volumes that are monitored to determine the radial distribution of each type of ion i, Ni(r) is large enough so that

Mixed-Salt Solutions of Polyanionic DNA

J. Phys. Chem. B, Vol. 103, No. 17, 1999 3493

TABLE 2: CMC Simulation Parameters for Solutions Containing a Model DNA Polyion (b ) 1.7 Å) and Mixed (2:1 and 1:1) Salts numbers of ionic species MC run number

C h 1:1 (mM)

11 12 13 14 15 15A 16 17 18 18A 19 20 21 22

16.5 16.5 16.5 16.5 72.0 72.0 72.0 72.0 72.0 72.0 239 239 239 239

cations univalent divalent 320 320 320 320 360 360 360 180 180 180 420 420 210 210

20 40 80 200 23 23 45 45 90 90 26 52 52 105

anions

polyion charges

approach distance a (Å)

cell radius Rc (Å)

cutoff radius Rd (Å)

fraction of accepted moves

240 280 360 600 286 286 330 210 300 300 412 464 284 390

120 120 120 120 120 120 120 60 60 60 60 60 30 30

10 10 10 10 10 13 10 10 10 13 10 10 10 10

229 229 229 229 118 118 118 118 118 118 100 100 100 100

150 130 95 80 65 65 55 50 50 50 55 50 45 45

0.572 0.582 0.611 0.566 0.571 0.556 0.566 0.656 0.616 0.570 0.577 0.559 0.644 0.593

its average value can be determined with an acceptably small statistical uncertainty ( 13 Å). Similar, although less pronounced, differences in the local extents of accumulation are observed for the univalent cation in the salt solutions where all of the cations are univalent (runs 6-10). Figure 3A indicates that C+MC(s) exceeds C+PB(s) by approximately 0.4 M, and Figure 3C shows that this difference decays to zero within the first 3-4 Å from the polyion surface. For either type of cation parts B and C of Figure 3 show that the excesses of C2+MC(r) and of C+MC(r) over the corresponding predictions of the PB equation are correlated with values of C-MC(r) that are significantly higher than those predicted by PB calculations, especially at the lower bulk salt concentrations (see Figure 1B). Ion-Polyion Preferential Interaction Coefficients for the Divalent Cation and the Univalent Anion in Model (DNAn-; M2+, X-) Solutions; Comparison with Results for the Univalent Cation. For model (DNAn-; M2+, X-) solutions, Γ2+

and Γ- were calculated independently by integrating the corresponding MC and PB radial ion distributions using eq 1d and 1a, respectively (cf. Figure 1). Results are tabulated in Table 4, and MC and PB predictions are compared in Figure 4. No experimental data are available for comparison. For both MC and PB calculations, the values of Γ2+ and Γ-, each independently determined by integration, are in all cases consistent with electroneutrality (2Γ2+ - Γ- ) 1). (Under all conditions investigated herein, ∑iziΓiMC deviates from 1 by less than 2% and ∑iziΓiPB by less than 0.5%.) Our primary focus in this paper is on values of Γi for the divalent and/or univalent cations, but values of Γi for the univalent anions, determined independently using the appropriate form of eq 1, also are presented. At sufficiently low salt concentrations, the signs and magnitudes of these molar scale preferential interaction coefficients are determined primarily by Coulombic interactions. (The contribution to Γi due solely to excluded volume, indicated explicitly in eq 1b, is relatively small under all conditions investigated herein.) At salt concentration(s) low enough so that, for example, any differences between molar scale and molal scale preferential interaction coefficients can be neglected, the theoretical range of Γ2+ in a (DNAn-; M2+, X-) solution is 0.33 < Γ2+ < 0.5, and the corresponding range of Γ- is -0.33 < Γ- < 0. For each type of ion, the lower bounds are the ideal values for a (hypothetical) electroneutral solution of a cylindrical polyanion with very low charge density in the presence of a uniform (random) distribution of M2+ and X- ions, and the upper bounds correspond to a (hypothetical) solution where divalent cation accumulation completely neutralizes the charge of the polyion and the presence of the polyion has no effect on the anion distribution.10 As predicted either by Manning’s counterion condensation (CC) theory27 or by the cylindrical PB equation,18 the low-salt limiting value of Γ- in a (DNAn-; M2+, X-) solution is -0.03, and therefore the corresponding limiting value of Γ2+ is 0.485. Figure 4 and Table 4 show that, when C2+bulk ) 10.4 mM, Γ2+MC and Γ-MC both are indistinguishable from their corresponding limiting values but detectably larger than Γ2+PB and Γ-PB, respectively. The limiting values of the Γi are free of inaccuracies introduced by the various kinds of approximations inherent in the analytical PB and CC formalisms,18,28 but they are not expected to remain applicable for a (DNAn-; M2+, X-) solution up to a bulk salt concentration of at least 10 mM, as our MC calculations indicate. This finding is particularly surprising because an analogous extended applicability of the limiting value of Γ+ is not predicted by our MC simulations on (DNAn-; M+, X-) solutions. For these systems the limiting law values of Γ+ and Γ- are 0.94 and -0.06, respectively. Table 4 and Figure 5 for the (DNAn-; M+,

Mixed-Salt Solutions of Polyanionic DNA

J. Phys. Chem. B, Vol. 103, No. 17, 1999 3497

Figure 3. (A) CMC and PB predictions of surface cation concentrations at the surface of model dsDNA as functions of bulk salt concentration for solutions containing only one type of salt. Values of C2+MC(s) (simulations 1-5) and of C+MC(s) (simulations 6-10) are plotted (open symbols) vs the logarithm of the bulk salt concentration and compared with the corresponding CiPB(s) (closed symbols). Differences between CMC and PB predictions of Ci(r) over the radial range are plotted where these differences are significant: (B) ∆C2+(r) (simulations 1 (0), 2 (]), 3 (n), 4 (g), 5 (4)); (C) ∆C+(r) (simulations 6 (]), 7 (g), 8 (4), 9 (O),10 (0)).

X-) solution show that these values significantly exceed both MC and PB predictions of Γ+ and Γ- at C1:1bulk ) 10 mM. With increasing C2:1bulk, Figure 4 shows that both Γ2+ and Γ- decrease (i.e., -Γ- increases). Each of these reductions is significantly smaller as calculated by MC than by PB, so that the difference ∆Γ2+ ≡ Γ2+MC - Γ2+PB increases from 0.023 ( 0.014 at C2:1bulk ) 10.4 mM to 0.048 ( 0.014 at C2:1bulk ) 223 mM (see Table 4). Discrepancies between Γ-MC and Γ-PB, which are twice as large as those between Γ2+MC and Γ2+PB, constitute a much larger percentage of Γ-MC, primarily because values of Γ- are much smaller in magnitude than are the corresponding Γ2+. Similar (but less pronounced) trends are exhibited by our calculations for (DNAn-; M+, X-) solutions (Figure 5 and Table 4). Both MC and PB values of Γ+ and Γdecrease with increasing 1:1 salt concentration, and the discrepancy between MC and PB predictions also may increase (slightly) with increasing salt concentration, from ∆Γ+ ) 0.016 ( 0.014 at C1:1bulk ) 10.6 mM to ∆Γ+ ) 0.031 ( 0.017 at C1:1bulk ) 210 mM.

Figures 4 and 5 illustrate an important but possibly counterintuitive aspect of thermodynamic counterion accumulation in these three-component (salt-polyion-solvent) solutions. For the high charge density cylindrical polyanion investigated herein, the cation preferential interaction coefficients Γ2+ and Γ+ are very large (near their respective maximum possible values of 0.5 and 1.0) at low bulk salt concentration and decrease with increasing salt concentration. Reasoning based on an analogy with chemical affinity (covalent or noncovalent) governed by mass action would lead one to expect these values to be small at low salt concentration and to increase with increasing salt concentration. Although Γ2+ and Γ+ are interpretable as the extents of cation accumulation per polyion monomer in the context of salt ion effects on processes involving a polyion,8 these coefficients are more fundamentally measures of the thermodynamic nonideality (nonrandomness) of the solution that is reflected by the nonuniformity of the local concentrations, as indicated by any of the forms of eq 1 (cf. Figures 1 and 2). For both di- and univalent cations, the steepness of each of these

3498 J. Phys. Chem. B, Vol. 103, No. 17, 1999

Ni et al.

Figure 4. CMC and PB predictions of Γi for 2:1 salt-DNA interactions. MC (0) and PB (9) results for (A) the divalent cation (Γ2+) and (B) the univalent anion (Γ-) are plotted as functions of C2+bulk.

Figure 5. CMC and PB predictions of Γi for 1:1 salt-DNA interactions. CMC (0) and PB (9) results for (A) the univalent cation (Γ+) and (B) the univalent anion (Γ-) are plotted as functions of C+bulk.

cation concentration gradients actually decreases with increasing bulk salt concentration, because the local cation concentrations near the polyion surface (Figure 3) are buffered against changes in bulk salt concentration, a phenomenon due to the screening of the polyion field by salt ions. The steepness of the anion concentration gradient also decreases with increasing salt concentration, even though the surface concentration of anions increases. Therefore, both cation and anion preferential interaction coefficients decrease in magnitude with increasing salt concentration. Although the physical exclusion of salt ions from the volume of solution occupied by the polyion also can affect the salt concentration dependence of the Γi, as indicated explicitly in eq 1b, such effects become significant only at salt concentrations higher than those investigated herein. The discrepancies between MC and PB values of cation preferential interaction coefficients for either (DNAn-; M2+, X- ) or (DNAn-; M+, X- ) solutions can be traced to the

Figure 6. CMC and PB radial distributions of ions for mixtures of 2:1 and 1:1 salts in a model dsDNA solution. Semilog plots in panels A-D (corresponding to results of MC simulations 15-18) display local (molar) concentrations of (A) divalent cations (C2+(r), 4), (B) univalent cations (C+(r), O), and (C) univalent anions (C-(r), 0) as functions of radial distance r from the axis of a cylindrical model dsDNA (b ) 1.7 Å, a ) 10 Å). Corresponding smooth curves are PB ion radial distributions calculated for each of the bulk cation concentrations determined in the CMC simulations.

differences between MC and PB predictions of local cations accumulated within a radial distance of 3 Å from the polyion surface (Figure 3, B and C). Integrations of these differences in local cation concentrations over this 3 Å range yield the values of ∆Ncationloc given in Table 4, which typically exceed the corresponding ∆Γcation by approximately a factor of 2. In the range a e r e a + 3, for cations CiMC(r) exceeds CiPB(r), but this trend is reversed in the next 5 Å and beyond so that the observed value of ∆Γcation is the composite of these two local effects. CMC and PB Radial Distributions of Ions in Mixed Cation (M2+, M+)-DNA Solutions. In Figure 6, the four panels (A-D) compare MC and PB radial distributions of M2+, M+, and X- ions for a series of divalent-univalent cation mixtures (DNAn-; M2+, M+, X- ). In each of these MC calculations (see Table 2), the average 1:1 salt concentration h 2:1MC was varied from 4.5 was fixed at C h 1:1MC ) 72 mM, and C mM (panel A; run 15) to 18 mM (panel D; run 18). For this

0.122 (0.003 0.112 (0.004 0.105 (0.002 0.099 (0.003 0.172 (0.003 0.164 (0.003 0.156 (0.002 0.154 (0.002 0.283 (0.003 0.280 (0.003 0.279 (0.002 0.277 (0.003 0.087 (0.011 0.065 (0.013 0.057 (0.009 0.045 (0.009 0.126 (0.010 0.103 (0.007 0.057 (0.009 0.035 (0.007 0.210 (0.011 0.192 (0.009 0.168 (0.007 0.148 (0.005 0.649 (0.003 0.430 (0.013 0.234 (0.004 0.140 (0.003 0.608 (0.017 0.440 (0.012 0.260 (0.007 0.158 (0.005 0.462 (0.006 0.331 (0.004 0.206 (0.002 0.121 (0.002 0.588 (0.013 0.349 (0.009 0.161 (0.009 0.073 (0.005 0.540 (0.012 0.349 (0.008 0.197 (0.009 0.119 (0.003 0.422 (0.018 0.312 (0.010 0.186 (0.007 0.101 (0.009 0.049 (0.018 0.066 (0.016 0.062 (0.016 0.058 (0.011 0.057 (0.020 0.077 (0.016 0.079 (0.013 0.083 (0.014 0.058 (0.016 0.058 (0.011 0.071 (0.013 0.087 (0.017 a

Reported range in ΓiPB results from uncertainty in corresponding MC Cibulk.

0.115 (0.010 0.229 (0.008 0.331 (0.005 0.401 (0.001 0.110 (0.009 0.198 (0.007 0.294 (0.003 0.345 (0.002 0.130 (0.005 0.197 (0.003 0.260 (0.002 0.303 (0.002 0.164 (0.008 0.295 (0.008 0.393 (0.011 0.443 (0.010 0.167 (0.011 0.275 (0.009 0.373 (0.010 0.428 (0.012 0.188 (0.010 0.255 (0.008 0.331 (0.011 0.390 (0.015 0.032 (0.010 0.267 (0.030 1.67 (0.015 7.33 (0.20 0.62 (0.10 2.47 (0.20 9.40 (0.38 25.4 (0.5 8.57 (0.47 20.5 (0.6 47.2 (1.4 105 (3.0 12.5 (0.2 14.1 (0.2 15.2 (0.2 15.7 (0.2 57.9 (0.6 62.3 (0.8 64.4 (0.7 66.2 (0.6 215 (2 218 (2 221 (2 224 (3 11 12 13 14 15 16 17 18 19 20 21 22

0.070 (0.010 0.073 (0.020 0.068 (0.013 0.071 (0.021 0.081 (0.010 0.078 (0.011 0.088 (0.012 0.086 (0.015 0.077 (0.016 0.078 (0.010 0.083 (0.015 0.089 (0.013

-Γ-PB a -Γ-MC Γ+PB a Γ+MC ∆Γ2+ Γ2+PB a Γ2+MC ∆N2+loc C2+bulk (mM) C+bulk (mM) MC run number

series, C+bulk was found to increase by CiPB(s), when comparisons are made at the same ion-polyion excluded volume. Unfortunately, no generally applicable, rigorous closedform relationship between excess surface ion concentrations and the distance of closest approach has been obtained. In principle, the appropriate value of a (or ao) can best be determined by analysis of the appropriate measurements on DNA-salt solutions. Fitting calculations of ΓiPB (ref 32) to the Donnan dialysis data of Strauss et al.33 for NaBr-NaDNA solutions (covering the range ∼10 mM to 1 M NaBr, but at relatively high DNA concentration) yields a distance of closest approach a ) 10 Å and, consequently, a radius ao = 8 Å if σ ) 2 Å. Using numerical solutions of the cylindrical PB equation to analyze small-angle neutron-scattering experiments in a mixed cation (tetramethylammonium, Na+) DNA solution34 yields a distance of closest approach a ) 14 Å, from which the authors estimate ao ) 8 Å. When the PB equation was used to analyze the effect of salt concentration on the Tm of DNA denaturation,32 equally good fittings were obtained for a range of choices of a for dsDNA provided that the difference ∆a ) ads - ass = 4 Å. These PB results, together with the minimum physically reasonable value implied by the primary structure of ssDNA, are consistent with the inference that ao is approximately 8 Å for dsDNA. In previous MC calculations, we used a cylinder radius of 10 Å for dsDNA and an ionic radius of σ ) 3 Å for univalent (e.g., Na+, K+, Cl-) ions.13,15,16,21 In the present study, the currently accepted values25,35 (ao ) 8 Å, σ ) 2 Å) are used for most of the calculations. To investigate the sensitivity of our results to different values of these structural parameters, some calculations were performed with the radii used previously. Table 6 compares MC and PB surface cation concentrations, local extents of cation accumulation (∆Ncationloc, the difference between MC and PB results for the number of cations per DNA phosphate within 3 Å of the DNA surface), and Γi for cations in representative solutions containing one or two salt components. Although CiPB(s) and CiMC(s) for both univalent and divalent cations are ∼25-50% smaller for a ) 13 Å than for a ) 10 Å, the corresponding Γi are j10% smaller for a ) 13 Å than for a ) 10 Å. Differences between MC and PB

0.328 ( 0.002

0.345 ( 0.002

0.385 ( 0.012

0.420 ( 0.012

0.075 ( 0.014

0.110 ( 0.009 0.167 ( 0.011

0.057 ( 0.014

0.087 ( 0.007 0.138 ( 0.010

0.057 ( 0.017

0.461 ( 0.001 0.469 ( 0.001 0.391 ( 0.002 0.420 ( 0.002 0.480 ( 0.011 0.492 ( 0.013 0.423 ( 0.013 0.459 ( 0.013

0.051 ( 0.020

predictions of Γi for cations are at most slightly smaller for a ) 13 Å than for a ) 10 Å. These comparisons are consistent with results reported by Sharp,36 who compared PB calculations of Γ with the MC calculations of Mills et al.13 for (DNAn-; M+, X- ) (a ) 10 Å, σ ) 3 Å). For various counterion radii σ and a constant cylinder radius ao ) 10 Å, Mills et al. found that the differences between CiPB(r) and CiMC(r) are greatly reduced when σ is increased from 0.5 to 3.5 Å. This observation was attributed to compensation between two opposing trends in effects on the higher order correlations that are neglected by the PB approximation built into eq 6: ion-ion Coulombic correlations and ion-ion excluded volume. Le Bret and Zimm26 reported increased differences between MC and PB counterion radial distributions calculated for salt-free DNA solutions where the counterions were assigned radii as large as 10 Å, which introduces a nonmonotonic structure into CiMC(r) near the DNA surface.

-0.039 ( 0.013 0.158 ( 0.005 0.119 ( 0.003

-0.033 ( 0.014 0.147 ( 0.004 0.114 ( 0.010

0.608 ( 0.017 0.540 ( 0.012

-0.068 ( 0.021

0.605 ( 0.012 0.538 ( 0.009

V. Comparisons with Previous Theoretical Calculations of Ion Distributions and Preferential Interaction Coefficients

0.086 ( 0.015

0.068 ( 0.011

0.081 ( 0.010

0.077 ( 0.011 0.082 ( 0.013 0.071 ( 0.010 0.083 ( 0.014 0.009 ( 0.008 0.036 ( 0.009 0.011 ( 0.008 0.040 ( 0.008 0.065 ( 0.010

a

Cbulk 2:1 ) 25.4 and Cbulk 1:1 ) 66.2 (18,18A)

Cbulk 2:1 ) 0.62 and Cbulk 1:1 ) 57.9 (15,15A)

Cbulk 1:1 ) 82.5 (8,8A)

Cbulk 1:1 ) 10.6 (6,6A)

Cbulk 2:1 ) 93.4 (3,3A)

Uncertainties determined as in Tables 4 and 5.

2.48 3.83 2.76 4.16 2.23 3.62 2.46 3.90 0.64 1.86 1.28 2.64 2.10 0.60 3.32 0.79 3.70 ( 0.13 5.22 ( 0.15 4.29 ( 0.14 5.81 ( 0.16 2.35 ( 0.07 3.96 ( 0.08 2.57 ( 0.08 4.47 ( 0.09 1.26 ( 0.09 1.78 ( 0.09 2.11 ( 0.04 2.51 ( 0.08 3.01 ( 0.13 0.55 ( 0.04 4.38 ( 0.18 0.70 ( 0.06 13 10 13 10 13 10 13 10 13 (2+) (+) 10 (2+) (+) 13 (2+) (+) 10 (2+) (+) ) 10.4 (1,1A) Cbulk 2:1

∆Nloc cation (ions/DNAP) a (Å)

surface cation concentration (M) MC PB bulk salt concentration(s) (mM) (run number)

TABLE 6: CMC-PB Comparisons for a ) 10 Å vs a ) 13 Åa

MC

0.869 ( 0.009 0.895 ( 0.010 0.757 ( 0.011 0.820 ( 0.013

Γ+

PB

0.861 ( 0.002 0.879 ( 0.002 0.748 ( 0.002 0.799 ( 0.002

∆Γ+

0.008 ( 0.011 0.016 ( 0.012 0.009 ( 0.013 0.021 ( 0.015

-0.067 ( 0.029

Γ2+

PB

0.019 ( 0.012 0.023 ( 0.014 0.032 ( 0.015 0.039 ( 0.015

Ni et al.

MC

∆Γ2+

3502 J. Phys. Chem. B, Vol. 103, No. 17, 1999

Ion Distributions. As reviewed recently by Jayaram and Beveridge,37 substantial effort has been devoted to comparisons between ion distributions predicted by the PB equation for rodlike polyions and by less approximate theories, such as those based upon the modified PB (MPB) equation, the hypernetted chain (HNC) equation, and/or canonical or grand canonical MC simulations. Calculations of ion-polyion distributions have been performed for a cylindrical model of dsDNA in both the presence and absence of added salt(s), where the cations are univalent and/or divalent. Physically appropriate ranges of cylindrical and ionic radii have been examined. In many calculations the DNA phosphate charges have been modeled by a continuous or discrete axial array, although a helical surface array, a grooved cylinder, and other modifications of the cylinder model have been considered. For a given model, properly equilibrated canonical or grand canonical MC simulations provide an exact description of the effects of Coulombic and excluded volume interactions on ionpolyion distributions (and hence, as discussed in the following subsection, on the thermodynamic properties that depend on ion-polyion interactions). For calculations of ion distributions (or the corresponding Γi), canonical MC simulations, like those performed herein, require significantly less computer time than do the corresponding grand canonical MC simulations for systems containing comparable numbers of particles. CMC simulations specify a priori the average concentration of each salt ion within the MC cell, but the corresponding bulk concentration must be obtained as an output subject to some level of statistical uncertainty. For dilute polyelectrolyte-salt solutions, comparisons of radial distributions, or thermodynamic properties, predicted by CMC and the various approximate analytical theories are performed at the same bulk salt concentration(s) in order to isolate the effects of salt ion interactions with a single polyion. For solutions containing one type of added salt (where the cation of the salt is the same as that introduced with the DNA) near the surface of a highly charged cylindrical polyion, the PB equation in general underestimates the local counterion (cation) concentration and overestimates the local anion concentration as compared to predictions based on the MPB or HNC equation or on MC simulations. The discrepancies are more severe for divalent than for univalent cations and for higher concentrations of added salt.26,35,38,39 Ion distributions calculated

Mixed-Salt Solutions of Polyanionic DNA by the MPB and HNC methods are almost equivalent to those calculated using MC simulations for the case of univalent cations and anions and a significant improvement over PB predictions for solutions containing divalent cations.26,35,38,39 Our MC-PB comparisons (Figures 1 and 2) are completely consistent with these observations. For DNA solutions, where both univalent and divalent cations are present (together with univalent anions), Rossky and co-workers38,39 calculated HNC and PB radial distributions to characterize the exchange of univalent and divalent cations near the model DNA over a range of bulk cation concentrations and examined the differences between HNC and PB distributions. General characteristics of our ion distributions (Figure 7) and our MC-PB comparisons for mixed cation cases are qualitatively consistent with the HNC-PB results of Rossky et al.38,39 and the MC-MPB studies of Das et al.35 A few studies have indicated that relatively minor effects on the Ci(r) result from introducing a helical array of polyion charge on the model cylinder26 or introducing grooves in the cylinder.25 For a cylinder radius of 8 Å, differences in Ci(r) calculated for the grooved cylinder and the impenetrable cylinder are confined to radial distances r e 12 Å, or 2 Å from the distance of closest approach.25 Thermodynamic Properties. Integrated ion radial distributions calculated from numerical solutions of the PB equation in various forms have been used to characterize the nonideality due to salt ion-polyion interactions in terms of various thermodynamic coefficients.18,20,31,36,40,41 In polyelectrolyte solutions containing only one type of salt component, equilibrium dialysis experiments can be used to evaluate the Donnan coefficient, which can be calculated on the basis of McMillanMayer theory as a kind of second-virial coefficient.19 This coefficient is related by standard thermodynamic transformations (as explained in the appendix to ref 9) to the corresponding salt-polyelectrolyte preferential interaction coefficient and hence to Γ-, calculated herein, using output from solutions of the PB equation (with eq 1a) or from CMC salt ion distributions (with eq 1d). Previously, we reported salt-DNA preferential interaction coefficients that were calculated using output from GCMC simulations on solutions containing only one type of excess salt (1:1). This approach requires a priori specification of the thermodynamic activity of the salt component (expressed as the appropriate product of ion densities).42 Both bulk and average ion concentrations are obtained from GCMC simulations as outputs subject to uncertainties due to statistical fluctuations. Only the average salt concentrations, in the presence and absence of the polyion at the same salt activity, are needed in order to calculate salt-polyelectrolyte interaction coefficients by the method introduced by Mills et al.13 Applying this method to dilute DNA solutions containing excess 1:1 salt, Mills et al.,13 Paulsen et al.,43 and Olmsted et al.16,44 have reported preferential interaction coefficients for a range of salt concentrations and structural model parameters. For conditions that overlap those investigated herein, we have confirmed that previously reported salt-DNA preferential interaction coefficients agree with the anion-polyion preferential interaction coefficients determined from our CMC anion radial distributions for DNA solutions containing 1:1 salt. Sharp36 compared univalent salt-DNA (co-ion) preferential interaction coefficients calculated for either the standard uniformly charged cylinder model or a structurally detailed “all atom” model. The results calculated using the 3D-PB equation for the latter model were systematically smaller (more negative by at least 5%) than those calculated using the 1D-PB equation for the cylindrical model. Surprisingly, this predicted discrep-

J. Phys. Chem. B, Vol. 103, No. 17, 1999 3503 ancy persists down to 1 mM, the lowest salt concentration investigated, even though the limiting (low-salt) value of the preferential interaction coefficient should be independent of the details built into the structural model of the polyion.18 Few previous calculations of preferential interaction coefficients or related thermodynamic quantities by any method have been reported for polyelectrolyte solutions containing divalent cations or divalent-univalent cation mixtures. For these systems, no previous MC calculations or directly comparable experimental measurements of Γi are available, although some published thermodynamic data for DNA solutions containing divalent cations can be interpreted using preferential interaction coefficients (Ni et al., in preparation) or related thermodynamic quantities.12 Using the cylindrical PB equation, Fogolari et al.11 investigated the distribution of divalent cations extending from the polyion surface to the bulk solution and, thereby, obtained a thermodynamic measure of divalent cation binding equivalent to our divalent cation preferential interaction coefficient (but lacking the excluded volume term in eq 1b). They used a Scatchard-type plot to analyze divalent cation binding as a function of divalent cation concentration and obtained Coulombic “cation binding” constants that exhibit the power dependence on C1:1 anticipated for cation-polyion interactions (cf. refs 1 and 8). Stigter and Dill12 used the PB potential due to the interactions of univalent cations and anions with DNA to predict the zero binding density extent of multivalent cation accumulation near DNA and, hence, obtained multivalent cation DNA binding constants without resorting to a Scatchard analysis like that presented by Fogolari et al.11 These binding constants also exhibit the anticipated power dependence on C1:1. Rouzina and Bloomfield45 have described multiply charged cation accumulation near DNA in solutions containing excess 1:1 salt on the basis of simple but approximate analytical expressions obtained from the cylindrical PB equation. Sharp36 and Chen and Honig46 reported PB electrostatic free energy calculations for mixed cation DNA solutions and for the binding of the λ repressor protein to DNA in a mixed cation solution, respectively. These PB calculations confirm the conclusions drawn from experiments on a variety of systems that even low concentrations of a divalent cation like Mg2+ greatly reduce the effects of changes in univalent salt concentration on oligocation and protein-DNA binding constants and on other thermodynamic properties of DNA solutions (see, for example, refs 2-4) and that the Debye-Hu¨ckel derived quantity “ionic strength” is inappropriate as a composition variable to describe these effects.47-49 Our present study quantifies the accumulation of divalent cations, such as Mg2+, using rigorous MC simulations on a simplified structural model of DNA, which previous studies (recently reviewed8) have indicated to be sufficient to capture quantitatively its polyelectrolyte effects on measurable thermodynamic properties. In a subsequent paper (Ni et al., in preparation), we will investigate the level of agreement between values of ΓiMC calculated for the cylindrical model of dsDNA and the experimentally observed salt concentration dependence of Mg2+-DNA interactions. Acknowledgment. We thank Dr. Arun Yethiraj for his advice on the CMC program, Dr. Jeff Bond for making available his code that solves the cylindrical PB equation in mixed counterion systems and Sheila Aiello for her assistance in preparing the manuscript. This research was supported by NIH Grant GM 34351 (M.T.R.). References and Notes (1) Record, M. T., Jr.; Zhang, W.; Anderson, C. F. AdV. Protein Chem. 1998, 51, 281.

3504 J. Phys. Chem. B, Vol. 103, No. 17, 1999 (2) Record, M. T., Jr.; deHaseth, P. L.; Lohman, T. M. Biochemistry 1977, 16, 4791. (3) Suh, W. C.; Leirmo, S.; Record, M. T., Jr. Biochemistry 1992, 31, 7815. (4) Capp, M. W.; Cayley, D. S.; Zhang, W.; Guttman, H. J.; Melcher, S. E.; Saecker, R. M.; Anderson, C. F.; Record, M. T., Jr. J. Mol. Biol. 1996, 258, 25. (5) Eisenberg, H. Biological Macromolecules and Polyelectrolytes in Solution; Clarendon: Oxford, 1976; p 29. (6) Schellman, J. A. Biophys. Chem. 1990, 37, 121-140. (7) Timasheff, S. N. Biochemistry 1992, 31, 9857. (8) Anderson, C. F.; Record, M. T., Jr. Annu. ReV. Phys. Chem. 1995, 46, 657. (9) Anderson, C. F.; Record, M. T., Jr. J. Phys. Chem. 1993, 97, 7116. (10) Record, M. T., Jr.; Richey, B. In ACS Sourcebook for Physical Chemistry Instructors; Lippincott, T., Ed.; American Chemical Society: Washington DC, 1988; p 145. (11) Fogolari, F.; Manzini, G.; Quadrifoglio, F. Biophys. Chem. 1992, 43, 213. (12) Stigter, D.; Dill, K. A. Biophys. J. 1996, 71, 2064. (13) Mills, P.; Anderson, C. F.; Record, M. T., Jr. J. Phys. Chem. 1986, 90, 6541. (14) Vlachy, V.; Haymet, A. D. J. J. Chem. Phys. 1986, 84, 5874. (15) Olmsted, M. C.; Anderson, C. F.; Record, M. T., Jr. Proc. Natl. Acad. Sci. U.S.A. 1989, 86, 7766. (16) Olmsted, M. C.; Bond, J. P.; Anderson, C. F.; Record, M. T., Jr. Biophys. J. 1995, 68, 634. (17) Allison, S. A. J. Phys. Chem. 1994, 98, 12091. (18) Anderson, C. F.; Record, M. T., Jr. In Structure and Dynamics: Nucleic Acids and Proteins; Clementi, E., Sarma, R., Eds.; Adenine Press: New York, 1983; p 301. (19) Schellman, J. A.; Stigter, D. Biopolymers 1977, 16, 1415. (20) Stigter, D. J. Colloid Interface Sci. 1975, 53, 296. (21) Mills, P.; Anderson, C. F.; Record, M. T., Jr. J. Phys. Chem. 1985, 89, 3984. (22) Allen, M. P.; Tildesley, D. J. Computer Simulation of Liquids; Clarendon: Oxford, 1987; p 156. (23) Fuoss, R. M.; Katchalsky, A.; Lifson, S. Proc. Natl. Acad. Sci. U.S.A. 1951, 37, 579. (24) Jayaram, B.; Beveridge, D. L. J. Phys. Chem. 1991, 95, 2506.

Ni et al. (25) Montoro, J. C. G.; Abascal, J. L. F. J. Chem. Phys. 1995, 103, 8273. (26) Le Bret, M.; Zimm, B. H. Biopolymers 1984, 23, 271. (27) Manning, G. J. Chem. Phys. 1969, 51, 924. (28) Fixman, M. J. Chem. Phys. 1979, 70, 4995. (29) Braunlin, W. H.; Strick, T. J.; Record, M. T., Jr. Biopolymers 1982, 21, 1301. (30) Paulsen, M. D. Ph.D. Thesis, University of WisconsinsMadison, 1998. (31) Klein, B. K.; Anderson, C. F.; Record, M. T., Jr. Biopolymers 1981, 20, 2263. (32) Bond, J. P.; Anderson, C. F.; Record, M. T., Jr. Biophys. J. 1994, 67, 825. (33) Strauss, U. P.; Helfgott, C.; Pink, H. J. Phys. Chem. 1967, 71, 2550. (34) Groot, L. C. A.; Kuil, M. E.; Leyte, J. C.; van der Maarel, J. R. C.; Cotton, J. P.; Jannink, G. J. Phys. Chem. 1994, 98, 10167. (35) Das, T.; Bratko, D.; Bhuiyan, L. B.; Outhwaite, C. W. J. Chem. Phys. 1997, 107, 9197. (36) Sharp, K. Biopolymers 1995, 36, 227. (37) Jayaram, B.; Beveridge, D. L. Annu. ReV. Biophys. Biomol. Struct. 1996, 25, 367. (38) Murthy, C. S.; Bacquet, R.; Rossky, P. J. J. Phys. Chem. 1985, 89, 701. (39) Bacquet, R.; Rossky, P. J. J. Phys. Chem. 1988, 92, 3604. (40) Marcus, R. A. J. Phys. Chem. 1955, 23, 1057. (41) Gross, L. M.; Strauss, U. P. In Chemical Physics of Ionic Solutions; Conway, B. E., Barradas, R. G., Eds.; Wiley: New York, 1966; p 361. (42) Torrie, G. M.; Valleau, J. P. J. Chem. Phys. 1980, 73, 5807. (43) Paulsen, M. D.; Richey, B.; Anderson, C. F.; Record, M. T., Jr. Chem. Phys. Lett. 1987, 139, 448. (44) Olmsted, M. C.; Anderson, C. F.; Record, M. T., Jr. Biopolymers 1991, 31, 1593. (45) Rouzina, I.; Bloomfield, V. A. Biophys. Chem. 1997, 64, 139. (46) Chen, S. W.; Honig, B. J. Phys. Chem. B 1997, 101, 9113. (47) Record, M. T., Jr.; Mazur, S. J.; Melancon, P.; Roe, J. H.; Shaner, S. L.; Unger, L. Annu. ReV. Biochem. 1981, 50, 997. (48) Record, M. T., Jr.; Ha, J. H.; Fisher, M. Methods Enzymol. 1991, 208, 291. (49) Lohman, T. M. Crit. ReV. Biochem. 1985, 96, 3.