Molecular Insights into the Ion-Specific Kinetics of Anionic Peptides

7098

J. Phys. Chem. B 2010, 114, 7098–7103

Molecular Insights into the Ion-Specific Kinetics of Anionic Peptides Joachim Dzubiella* Physics Department T37, Technical UniVersity Munich, 85748 Garching, Germany ReceiVed: February 4, 2010

The action of NaCl vs KCl on the static and kinetic behavior of a fully charged and unfolded polyglutamic acid (PGA) chain is investigated by extensive explicit-water computer simulations. Ion-specific shrinking of the PGA coil size with increasing salt concentration is observed and is consistent with intrinsic viscosity measurements. The PGA relaxation kinetics is found to be nearly exponential in KCl on a Rouse/Zimm time scale (=1 ns), whereas NaCl induces a 10- to 100-times slower, highly nonexponential relaxation. The slow decay can be traced back to Na+ ions bridging anionic groups with scale-free power-law residence time distributions. This “transient cross-linking” may explain cation-specific slowing down of (bio)polymer kinetics observed in a variety of experiments. A systematic test using different force-field combinations in the simulations corroborates the qualitative trends, while quantitatively, the kinetic rates in the NaCl simulations significantly depend on the particular choice of water and ion parameters. I. Introduction Protein folding and assembly kinetics are typically interpreted in terms of a diffusional search through reduced energy landscapes.1-5 The total effective diffusion (or friction) profile along a suitable reaction coordinate arises from both peptidesolvent and intrapeptide interactions, suggesting a conceptual separation between solvent-induced and internal friction contributions, respectively. Recent experiments have demonstrated that internal friction constitutes a major contribution to the total friction, and solvent viscosity is not the only quantity that governs protein kinetics.6,7 However, solvent-induced and internal friction mechanisms are intimately connected; with addition of guanidinium ions, for instance, intrachain friction in unfolded glycine-serine peptide chains decreases, although solvent viscosity increases.8 This “lubrication” effect has been attributed to weak binding of guanidinium to the peptide backbone, thereby breaking intrapeptide hydrogen bonds and allowing the peptide parts to slip off each other. The investigation of the ion-specific action on proteins (Hofmeister effects) has a long history,9 but the molecular mechanisms are still under exploration. In a recent series of experiments, it has been shown that even the simple cations sodium (Na+) and potassium (K+) exhibit considerably different behavior in the interaction with protein surfaces, where Na+ is favored over K+.10-12 One apparent reason is the stronger attraction of sodium to acidic (anionic) surface groups, in particular carboxylate groups. Although these static properties have received much attention, not much is known about their consequences for biomolecular kinetics. Experimental hints have been given in studies of Na+- and K+-specific polyglutamic acid aggregation kinetics,13 folding kinetics of halophilic (“saltloving” and very acidic) proteins,14,15 or DNA.16 Similar specific relaxation has been observed in polymer melts.17,18 In contrast to the guanidinium cation,8 a slowing down of the kinetics for Na+ was always observed. A molecular understanding of ionspecific (bio)molecular kinetics is still lacking. In this study, we investigate the cation-specific (Na+ vs K+) relaxation kinetics of a short and highly charged polyglutamic * E-mail: [email protected].

acid (PGA) chain by all-atom computer simulations including explicit water molecules. On one hand, PGA may serve as a model peptide for partly acidic protein fragments or other strongly cation-binding polyelectrolytes. On the other hand, sequences of PGA occur in the acidic tails of transcription factors19 and have numerous industrial applications, such as wastewater treatment, food products, cosmetics, drug delivery, medical adhesives, vaccines, and tissue engineering20 and thus are of practical interest in their own right. The ion-specific effects of Na+ and K+ on (macroscopic) PGA properties are well documented and include electrical conductivity,21 ion permeability and swelling in cross-linked membranes,22 sedimentation velocity,23 helix-coil transition,24 and intrinsic viscosity measurements.25 In this contribution, we rationalize and provide molecular insight why the kinetic and frictional properties of PGA chains (and related polyelectrolytes) specifically change upon the addition of simple cations. Our findings may help to interpret the ion-specific kinetics found in experiments,13,16-18 in particular, for halophilic proteins.14,15 II. Methods Our molecular dynamics (MD) computer simulations are performed using the parallel module pmemd in the simulation package Amber9.0 with the ff03 force field for the peptide and the rigid and nonpolarizable TIP3P water model for the solvent.26,27 We use the NpT ensemble with N = 6000 atoms, pressure p ) 1 bar, and temperature T ) 300 K by coupling to a Berendsen barostat and Langevin thermostat,26 respectively. The cubic and periodically repeated simulation box has edge lengths L = 40 Å, including ∼2000 water molecules. Electrostatic interactions are calculated by particle mesh Ewald summation, and all real-space interactions have a cutoff of 9 Å. The peptide is generated and its trajectories are analyzed using the tleap and ptraj tools in the Amber package,26 respectively. All simulation snapshots are visualized using VMD.28 We investigate a 12 amino acid long PGA chain with the acetyl (Ace) and amine (Nme)-capped sequence Ace-(Glu)12Nme. We simulate at the fully deprotonated state, mimicking a neutral pH, so that the peptide has a net charge of -12e. The

10.1021/jp1010814  2010 American Chemical Society Published on Web 05/03/2010

Ion-Specific Kinetics of Anionic Peptides

J. Phys. Chem. B, Vol. 114, No. 20, 2010 7099

TABLE 1: Atom-Atom LJ Parameters and Partial Charges, q, Used in Our Worka ion

σ/nm


/(kJ/mol)

Na Na+ K+ Na+ (SPC/E) Na+ (TIP3P) ClCl- (SPC/E) Cl- (TIP3P) SPC/E O H TIP3P O H

0.2584 0.2350 0.3320 0.2160 0.2439 0.4401 0.4830 0.4477

0.4184 0.5440 0.4185 1.4755 0.3659 0.4184 0.0535 0.1489

+1 +1 +1 +1 +1 -1 -1 -1

0.3169

0.6501

-0.8476 +0.4238

0.3151

0.6364

-0.8340 +0.4170

+

charge q/e

ref 30 35 30 38 38 30 38 38 34 27

a The ions and water model set in bold font define the primary investigated force field.

ions Na+, K+, and Cl- are modeled as nonpolarizable LennardJones spheres. The charge and interaction parameters for our primary force field are provided by Dang30 (see Table 1). We simulate salt concentrations ranging from c ) 0.08-4 M and, as reference simulations, the counterion-only (CO) case using 12 Na+ ions and a simple neutralizing background (NB) with no explicit ions at all. For c = 2 M, we have tested the influence of five other force-field combinations as described in the Appendix. We equilibrate the systems for =50 ns before gathering statistics in an =0.5-µs-long simulation for each particular system. The PGA chain freely moves in the simulation box. To check for finite size effects, which could be anticipated for no or very weak screening by salt, we simulated selected systems in a bigger box with edge lengths L = 50 Å. The results are briefly discussed in the Appendix. We do not observe any appreciable folding into secondary structures induced by the salt in accordance with experiments of fully ionized PGA chains.24 We find only very rare folding events into partial R-helical states, whereas the average helicity stays below 5% in most cases. Only for 0.32 M NaCl was the helicity enhanced to 13% on average. In a recent simulation study, a higher helicity was measured at a similar NaCl concentration for a somewhat longer PGA chain.29 III. Results A. Coil Extension. A simulation snapshot of the pure PGA chain is shown in Figure 1 a). The PGA’s end-to-end distance bNme - b RAce| is defined to be the spatial separation between R ) |R the centers of mass of the Ace and Nme caps of the peptide. As a reference system, we have simulated the PGA chain without any explicit counterions, just using a neutralizing background (NB). The end-to-end distance distribution P(R) for this system is plotted in Figure 2 and exhibits a shifted Gausslike behavior. This is expected for a short, unfolded (semiflexible) polymer chain31 that is not collapsed due to charge repulsion between the side chains. To estimate the width, σ, of the distribution, we apply a skewed and shifted Gaussian fit8 of the form P(R) ) aR2 exp[-(R - b)2/σ2], with fitting constants a and b. We obtain σ ) 5.7 Å for the PGA chain in the NB. As a measure of PGA chain size, we calculate the mean end-toj ) (〈R2〉)1/2. The persistence length, end distance defined by R j 2 ) 2lclp - 2lp2[1 - exp(-lc/lp)], lp, can be extracted via31 R where lc ) 45.5 Å is the contour length of our considered PGA j ) 19.7 Å and lp = 5 Å, the latter chain. For the NB, we find R indeed in the range of unfolded PGA and other peptide chains.32

Figure 1. MD snapshots of the PGA chain. (a) Extended chain configuration. Anionic groups are rendered in red. Water molecules and ions are omitted for clarity. (b) Long-lived collapsed chain configuration in NaCl. The blue spheres represent Na+ ions, which are dynamically selected such that they have been bound to the peptide for more than 1 ns.

Figure 2. Probability distribution P(R) of the end-to-end distance R of the PGA chain in different salt solutions. MD results (symbols) are fitted with a skewed and shifted Gaussian (lines) (see text).

As an alternative measure for the chain size, we evaluate the radius of gyration, Rg, and we find Rg ) 10.0 Å in the NB. Note, however, that in the NB case, electrostatic interactions between the PGA images are not screened in the simulation, and the chain size depends on the box size (see the system size discussion in the Appendix). The effects of the addition of counterions and salt on the P(R) distribution are presented in Figure 2. We find Gaussianj, like distributions for all investigated systems for which σ, R and Rg are summarized in Table 2. Upon insertion of neutralizing j and Rg decrease by roughly 10%. Further Na+ counterions, R addition of NaCl amplifies the decrease in the chain size. At a j and Rg are shrunk by roughly 22%, concentration of =2 M, R and the width, σ, is half that of the PGA chain in the NB. The decrease in the PGA chain size with salt concentration is considerably smaller if Na+ is replaced by K+. The addition of

7100

J. Phys. Chem. B, Vol. 114, No. 20, 2010

Dzubiella

TABLE 2: Properties of the PGA Chain in Neutralizing Background (NB), Counterions Only (CO) and Different Salt Concentrations, ca NB CO NaCl

KCl

c/M

σ/Å

j /Å R

Rg/Å

nc (COO-)

nc (CO)

τ/ns

0.00 0.00 0.08 0.16 0.32 1.98 3.97 0.3 1.92 3.76

5.7 5.7 5.7 5.3 4.2 3.2 5.6 5.7 5.7 5.5

19.7 18.4 17.5 16.9 16.2 15.1 16.1 19.0 18.4 18.4

10.0 9.0 8.7 8.6 7.4 7.8 8.1 9.9 9.3 9.0

0.46 0.54 0.56 0.72 1.24 1.70 0.24 0.48 0.70

0.13 0.20 0.16 0.28 0.28 0.42 0.07 0.11 0.18

0.6 4.4 6 7 10 11 9 0.4 0.5 0.6

a The variable σ is the width of the end-to-end distance j is the mean end-to-end distance, Rg is the radius distribution P(R), R of gyration, nc is the average cation coordination number around COO- and CO groups, and τ is the typical relaxation time of the end-to-end distance autocorrelation function.

KCl decreases the chain size by only 7% at 1.9 M and by 10% at 3.8 M relative to the neutral NB case. Thus, whereas shrinking of the chain extensions could have been anticipated due to simple electrostatic screening by salt, the magnitude of the effect is strongly dependent on the type of cation. The ion-specific compression of the PGA chain is in qualitative agreement with experiments on the intrinsic viscosity, η, of fully ionized PGA coils.25 From polymer theory, it is well established that the intrinsic viscosity scales with chain size as η ∝ Rg3.31 In agreement with the enhanced shrinking of PGA coils in NaCl vs KCl solutions found above, the experiments measure a larger decrease in η for NaCl.25 For a rough quantitative comparison, in the molar range 1-2 M, the experiments yield an ∼20% larger decrease in η in NaCl vs KCl, corresponding to a 6% smaller Rg of the coil in NaCl. In our simulations, the PGA coil compression in NaCl seems somewhat overexpressed as we find a decrease of Rg by 16%. Note, however, that the experiments are performed in the limit lc . lp and high dilution of the PGA chains, in contrast to our simulations. In addition, we will show later that the exact value of the PGA chain size depends also on the particular choice of the NaCl force field. The qualitative larger decrease in Rg for NaCl vs KCl, however, is robust for all investigated force fields. B. Ion Binding. The ion-specific compression of the PGA chain was loosely attributed to the individual binding of the cations to anionic peptide groups.25 To investigate this notion in more detail, we have calculated radial distribution functions (RDFs) between the cations and the electronegative peptide atoms, that is, the backbone oxygen (Obb) and side chain oxygen (Osc). Examples for a salt concentration of =2 M are shown in Figure 3. Clearly, Na+ binds much more strongly to the oxygen atoms than to K+. In particular, the interaction between the carboxylate Osc and Na+ is considerable, expressed by a RDF contact peak of more than 20. This trend of a higher Na+ affinity vs K+ to carboxylates and carbonyls is in agreement with recent experimental and computational studies on ion specificity at protein surfaces.10-12 The average coordination number, nc, of the cations in the first hydration shell of the oxygen atoms is also given in Table 2. The coordination of Na+ is always a factor 2-3 higher than for K+, in quantitative agreement with cation binding to protein surfaces.10 The observation of specific binding raises the question, what are the typical residence times of the bound cations? The residence time distributions Pb(t) are shown in Figure 3a for both cations binding to the carboxylate oxygen. Binding is

Figure 3. The RDF between Na+ and K+ and the side chain carboxylate oxygen (Osc) or the backbone carbonyl oxygen (Obb) at c = 2 M. (Inset a) Log-log plot of the residence time distribution, Pb(t), of Na+ and K+ in the first solvation shell of Osc. K+ binding decays exponentially with a time constant τ = 45 ps; the Na+ ion exhibits a power law behavior Pb(t) ∝ t-2. (Inset b) The RDF between the headgroups of glutamic acids (Glu) side chains in NaCl and KCl solutions at a salt concentration of =0.3 M.

defined by a mutual distance closer than the first minimum in the cation-oxygen RDF, which is 3.15 and 3.60 Å for Na+ and K+, respectively. All K+ distributions exhibit an exponential decay Pb(t) ∝ exp(-t/τ), and we find a typical binding time of τ = 45 ps. In striking contrast, the Na+ residence time distribution obeys a power law behavior Pb(t) ∝ tR for long times (t J 0.1 ns) with an exponent R ) -2 (see Figure 3a). We observe the same behavior and exponent for the binding of cations to the backbone carbonyls (not shown). The power-law behavior points to a time-scale free binding of sodium, and long residence times in the nanosecond regime are possible. A power law residence time distribution with an exponent of R ) -1.6 has also been detected for Na+ binding to lipid carbonyls in a simulation study of a lipid bilayer system.33 The power law behavior in our study is found to be very robust against changes in the force field. Why does the cationic binding lead to the compression of anionic polyelectrolytes in the first place? Inspection of our MD trajectories shows that the strong Na+ binding leads to cationinduced bridging of two or more anionic peptide groups. In other words, one or two bound Na+ ions can be shared between carboxylic groups what eventually leads to their effective attraction. As a rough quantification, we have calculated the RDF between the charged headgroups of the (Glu) side chains (see the inset b in Figure 3). Whereas for all KCl concentrations, the RDF does not show significant structure above unity, that is, no mutual attraction, the RDF in NaCl exhibits a contact peak of about 3, indicating a significant effective attraction on the order of the thermal energy kBT ln(3) = kBT. We observe the same effect, but less strong, between carboxylates and backbone carbonyls (not shown). An exemplifying MD simulation snapshot of Na+ ions bridging peptide fragments is shown in Figure 1b, where the PGA chain forms a looplike structure around the cations. In view of the time-scale free binding distribution of Na+, these chain configurations can be assumed on a long, nanosecond time scale. This temporary cross-linking effect explains the ion-specific chain compression, and an increased internal peptide friction must be anticipated. These cross-links may also have structural relevance in biomolecular

Ion-Specific Kinetics of Anionic Peptides

J. Phys. Chem. B, Vol. 114, No. 20, 2010 7101

Figure 4. End-to-end distance autocorrelation function CRR(t) for the PGA peptide vs time, t, in a log-linear plot. KCl yields an exponential decay with τ = 0.5 ns, and NaCl gives a nonexponential decay that is at least 1 order of magnitude slower (see Table 2 for numbers). Inset: Survival probability S(t) of noncontact states of the two end monomers (Ace and Nme) in =2 M KCl and NaCl.

stabilization because glutamic acid residues locked by bound Na+ ions have been found in crystal structures of halophilic proteins.14 C. PGA Kinetics. To exemplify the ion-specific relaxation kinetics of the PGA chain, we focus on the first relaxation mode; that is, the kinetics of of the end-to-end distance R. We define the autocorrelation of the latter in time t to be

CRR(t) ) 〈R′(t) R′(0)〉/〈R′(0)2〉

(1)

j is shifted by R j so that the function where R′(t) ) R(t) - R CRR(t) decays to zero for convenience. The autocorrelation function for selected systems is presented in Figure 4. The reference system without ions and the KCl systems exhibit a decay well described by a single exponential CRR(t) ∝ exp(-t/ τ) resembling Rouse-like dynamics of ideal polymers.31 We calculate the relaxation time scale by integration τ ) ∫0∞ CRR(t) dt and find τ = 0.5 ( 0.1 ns for the chain in KCl. This is on j 2/ the same order of magnitude as the Rouse time τR ) NbR 2 3π D = 1 ns, where we assumed Nb ) 13 bonds and a typical molecular diffusion constant on the order of D = 10-5 cm2/s. We find that the same analysis for a Zimm chain31 yields the same order of magnitude of ∼ 1 ns. From that perspective, the PGA chain in KCl solution resembles a simple ideal polymer. In strong contrast, the relaxation behavior of the PGA chain in the Na+ environment differs significantly from KCl: already in the counterion-only (CO) case, a nonexponential relaxation is visible for long times t J 0.5 ns. The typical relaxation time from integration over the distribution is τ = 4.4 ns, 1 order of magnitude larger than for the KCl solutions. Increasing NaCl concentration amplifies this effect, as is visible in Figure 4. The relaxation times are all summarized in Table 2 and range from 6 to 11 ns, which are a factor of roughly 20 larger than in the KCl solution. The PGA in the NaCl solution thus exhibits a behavior qualitatively different from a Rouse/Zimm chain. Note that ideal behavior predicts decreasing relaxation times for smaller chain extensions,31 opposite to the observed behavior. Another experimentally accessible kinetic quantity is the endto-end loop formation rate.8 This rate can be estimated from the survival probability, S(t), of noncontact states of the two end monomers, Ace and Nme. We define them to be in noncontact if R > 6.5 Å, which constitutes a minimum in the

probability distributions P(R) (see Figure 1). We assume that every saved noncontact configuration R > 6.5 Å for each MD trajectory run is a valid starting point in the analysis. The survival probability S(t) averaged over all starting points is plotted in the inset to Figure 4 for the NB case and =2 M KCl and NaCl. The loop formation rates, k-1 ) ∫0∞ S(t) dt, are 1/(24 ns), 1/(14 ns), and 1/(109 ns) for NB, KCl, and NaCl, respectively. Consistent with the slow segmental relaxation in NaCl, the loop formation kinetics is about 1 order magnitude slower in NaCl than in KCl. On the basis of our analysis, the increased friction and the new time scales must be attributed to temporary cross-linked (trapped) chain fragments induced by the Na+ ions with a power law residence time distribution. It is, indeed, well established that systems with multiple trapping or other manifestations of disorder can lead to anomalous kinetics.41 In the investigated strongly anionic peptide, the appearance and magnitude of trapping is controlled by the nature of the cation. The observed “stretched” exponentials resemble the slow relaxation in glassforming liquids.42 At a critical concentration of a certain type of salt, an ion-induced kinetic arrest, that is, glassification or gelation, of internal peptide kinetics may be triggered. We actually find strong indications of glassy behavior for other NaCl force fields, which induce even stickier interactions (see the Appendix). IV. Summary and Concluding Remarks In summary, we have given molecular insight into cationspecific effects on the kinetics of an anionic, unfolded peptide. The previously observed stronger binding of Na+ vs K+ ions to electronegative peptide groups10-12 induces effective attractions between peptide fragments by temporary bridging and leads to cation-specific chain compression consistent with experiments.25 The kinetic consequence is an increased internal friction leading to a slow, nonexponential relaxation of the dynamic modes and decreased loop formation rates in NaCl solution. Our microscopic insight may help interpret experimentally observed cation-specific slowing down of (bio)polyelectrolye kinetics,13,16 in particular, for halophilic proteins.14,15 Similar mechanisms may be at work in polymer melts.17,18 Qualitative trends survive for a variety of force-field combinations, but the kinetic rate numbers in the NaCl simulations significantly depend on the particular models employed. The slow cosolute-induced segment relaxation and its implications for protein (folding) dynamics promises to be highly interesting for future considerations. On a more coarse-grained level, recent work has demonstrated that some sort of attractive interactions between polymer monomers yields anomalous kinetics of the relaxation modes43 and may dominate the dissipated work in protein unfolding.44 In general, anomalous diffusion may be a common feature in protein and peptide dynamics43,46 and challenges standard (normal diffusion) Smoluchowski approaches to protein folding.45 Our work shows that anomalous diffusion can also be induced and controlled by cosolute binding. The mechanisms behind the change of protein kinetics, however, may qualitatively depend on the particular and specific nature of the interaction. The guanidinium cation, for instance, has been shown to decrease friction in neutral (GlySer)n peptides.8 Apparently, (GlySer)n in pure water is in a collapsed state with high internal friction exhibiting anomalous kinetics.46 The “lubrication” mode of action of guanidinium has been attributed to the breaking of peptide-peptide hydrogen bonds by weak cosolute binding, which apparently leads to normal diffusion.

7102

J. Phys. Chem. B, Vol. 114, No. 20, 2010

Dzubiella

TABLE 3: Same as Table 2 for Other Investigated Force Fields and Systems Sizesa test case 1 (NaCl) 2 (NaCl) 3 (NaCl) 4 (NaCl) 5 (KCl) NaCl primary KCl primary NB NB (big) CO (big) NaCl (big) NaCl (big) a

force field

c/M

σ/Å

+ 30

1.98 1.98 1.98 1.98 1.92 1.98 1.92 0.00 0.00 0.00 0.07 1.87

4.0

(Na ) + SPC/E (Na+)35 + SPC/E (Na+)38 + TIP3P (Na+)38 + SPC/E (K+)30 + SPC/E (Na+)30 + TIP3P (K+)30 + TIP3P SPC/E TIP3P (Na+)30 + TIP3P (Na+)30 + TIP3P (Na+)30 + TIP3P

2.1 5.7 3.2 5.7 5.7 5.7 5.7 5.7 3.4

j /Å R

Rg/Å

nc (COO-)

nc (CO)

τ/ns

16.9 12.9 15.2 13.8 19.1 15.1 18.4 20.3 25.4 20.6 17.9 15.6

8.1 8.5 7.4 8.1 9.8 7.8 9.3 11.1 10.8 9.3 8.8 7.6

0.59 0.71 0.65 0.78 0.33 1.24 0.48

0.32 0.47 0.22 0.42 0.09 0.28 0.11

0.35 0.43 1.20

0.11 0.15 0.30

4.8 J100 30 J60 0.5 11 0.5 0.5 0.4 4.0 5.5 10

The primary investigated force fields for c = 2 M are given (in bold font) for direct comparison.

In our study, the opposite mechanism is uncovered, in which binding by “hard” ions increases internal friction by transient cross-links and leads to anomalous kinetics. We are hopeful that our work sheds more light onto cosolute binding effects on protein folding and assembly kinetics. More experimental studies, for example, by time-resolved Fo¨rster resonance energy transfer (FRET), are highly desirable.8 In particular, the action of complex denaturants such as guanidinium and urea deserve further attention, and systematic studies on specific salt effects should follow. We expect a similar strong influence of other strongly binding cations on anionic peptides, such as lithium, or polyvalent cations, such as Mg2+ or Ca2+. Large effects may also be anticipated by exchanging the anion that has been found to considerably alter the unfolding kinetics of a halophilic protein.14,15 On the basis of our results, one may expect that a possible anion effect is to influence the strength of the cation binding to the peptide. On the other hand, a weak anion, such as iodide or thiocyanate, may directly compete with the cation in changing internal friction due to lubrication. Acknowledgment. J.D. is grateful to C. de Oliveira, A. Erbas, Y. von Hansen, M. Hinczewski, R. R. Netz, and B. Rotenberg for useful discussions; the Deutsche Forschungsgemeinschaft (DFG) for support within the Emmy-Noether-Program; and the Leibniz Rechenzentrum (LRZ) Mu¨nchen for computing time on HLRB II. Appendix Force Field and System Size Dependence. The phenomena discussed in our work result from strong and specific binding effects of cations on the peptide carbonyl and carboxylate oxygen atoms. One can easily argue that these effects might be sensitive to small changes in the underlying (classical) force fields that are used in the MD simulation, especially for the stronger binder Na+. Therefore, we have tested four other force fields for the NaCl solution and one other force field for KCl at a concentration c ) 1.98 M, in which different parameters {ε, σ} for the Lennard-Jones interaction

[( σr ) - ( σr ) ]

VLJ ) 4ε

12

6

(2)

and a different water model is used. All employed parameters are summarized in Table 1. The primary parameters (NaCl and KCl from Dang30 and TIP3P water, which is the default water model in AMBER) used in the main part of this work are shown in bold font. In test case 1, we have replaced TIP3P water27 in the NaCl solution with SPC/E water,34 as in Dang’s work.30 In

test case 2, we have replaced Na+ from ref 30 with another nonpolarizable Na+ from ref 35 that has been parametrized in SPC/E water. This force field has been recently used to model monosodium glutamate solutions, and results consistent with experiments have been found for glutamate structure and hydration numbers.36 Both test cases 1 and 2 are also very similar to very recently proposed force fields for sodium in SPC/E water based on thermodynamic solvation properties (see sets 5a and 5b in ref 37). In test cases 3 and 4, we have used other very recently developed NaCl force fields for biomolecular simulations in SPC/E and TIP3P water,38 respectively. In test case 5, we use Dang’s K+ in SPC/E water instead of TIP3P. We have not tried to use the Amber9.0 default parameter for NaCl and KCl; they exhibit unreasonable ion clustering due to the combination of Åqvist’s cations39 and Dang’s anions, as has been recognized in the literature.38,40 All employed force fields in this work use the Lorentz-Berthelot mixing rules, in which σij ) (σii + σjj)/2 and εij ) (εiiεjj)1/2 for atoms i and j. As a reference, we have also simulated the NB system in SPC/E water. For all five test systems, we have performed the same analysis as for the primary force field. The results are summarized in Table 3. We find that the qualitative trends remain the same, independent of the force field. The mean end-to-end distance, j , and the radius of gyration, Rg, are smaller in NaCl than for R KCl at the same concentration, indicating cation-specific compression of the PGA chain. Quantitatively, the radius of gyration is between 20% (test case 3) and 9% (test case 2) smaller than that of the PGA chain in KCl. For all test cases, a much larger relaxation time scale τ is also observable for NaCl vs KCl, showing that the cation-specific trend in the peptide kinetics is robust against small changes in the force field. Quantitatively, however, the value of τ is sensitive to the force field and ranges from τ ) 4.8 ns (test case 1), which is still 1 order of magnitude slower than for KCl, to more than 60 and 100 ns for test cases 4 and 2, respectively. These long relaxation times comparable to the total simulation time may indicate a kinetic arrest, and τ may even grow for longer simulation times. This “glassy” dynamics in cases 2 and 4 must be attributed to an even stronger binding to the carboxylate and carbonyl groups than for the other force fields, as reflected in the average coordination numbers, nc (see Table 3), and stickier interactions between the peptide groups are induced. In addition, for test cases 2 and 4, no estimate of the width σ of the P(R) distribution was possible because of a highly non-Gaussian form of P(R). Kinetics traps in the sampling due to the high stickiness of the interactions may be responsible for this. The power law residence time distribution Pb(t) ∝ t-2

Ion-Specific Kinetics of Anionic Peptides of Na+ to peptide oxygens can be found for all tested force fields. All trends are qualitatively the same and consistent with experiments, and thus corroborate our main conclusions. Quantitative predictions of ionic force fields in biomolecular simulations, however, should be considered critically. We hope that which force field of the above comes closest to reality will be illuminated in not-too-far future studies in combined compuational and experimental efforts. Regarding finite size effects, we find a larger extension of the highly charged chain in a bigger box in the case of no screening or very weak screening by counterions only (CO) (see Table 3). For the lowest salt concentration considered, c ) 0.07 M; however, the screening length becomes as small as =10 Å, and the results slowly converge to the one of the smaller box size. For c ) 1.87 M, we indeed find very similar behavior for both box sizes considered. Therefore, the trends with salt type and concentration do not change with the computational box size. References and Notes (1) Bryngelson, J.; Wolynes, P. J. Phys. Chem. 1989, 93, 6902. (2) Camacho, C.; Thirumalai, D. Proc. Natl. Acad. Sci. 1993, 90, 6369. (3) Socci, N.; Onuchic, J.; Wolynes, P. J. Chem. Phys. 1996, 104, 5860. (4) Dill, K.; Chan, H. Nat. Struct. Mol. Biol. 1997, 4, 10. (5) Hummer, G.; Garcia, A.; Garde, S. Phys. ReV. Lett. 2000, 85, 2637. (6) Pabit, S.; Roder, H.; Hagen, S. Biochemistry 2004, 43, 12532. (7) Cellmer, T.; Henry, E.; Hofrichter, J.; Eaton, W. Proc. Natl. Acad. Sci. 2008, 105, 18320. (8) Mo¨glich, A.; Joder, K.; Kiefhaber, T. Proc. Natl. Acad. Sci. 2006, 103, 12394. (9) Baldwin, R. L. Biophys. J. 1996, 71, 2056. (10) Vrbka, L.; Vondra´sek, J.; Jagoda-Cwiklik, B. Proc. Natl. Acad. Sci. 2006, 103, 15440. (11) Uejio, J.; Schwartz, C.; Duffin, A.; Drisdell, W.; Cohen, R.; Saykally, R. Proc. Natl. Acad. Sci. 2008, 105, 6809. (12) Aziz, E.; Ottosson, N.; Eisebitt, S.; Eberhardt, W.; Jagoda-Cwiklik, B.; Va´cha, R.; Jungwirth, P.; Winter, B. J. Phys. Chem. B 2008, 112, 12567. (13) Colaco, M.; Park, J.; Blanch, H. Biophys. Chem. 2008, 136, 74. (14) Madern, D.; Zaccai, G. Extremophiles 2000, 4, 91. (15) Bandyopadhyay, A. K.; Krishnamoorthy, G.; Padhy, L. C.; Sonawat, H. M. Extrempophiles 2007, 11, 615. (16) Gray, R. D.; Chaires, J. B. Nucleic Acids Res. 2008, 36, 4191. (17) Bergman, R.; et al. Electrochim. Acta 1995, 40, 2049.

J. Phys. Chem. B, Vol. 114, No. 20, 2010 7103 (18) Mos, B.; Verkerk, P.; Pouget, S.; vanZon, A.; Bel, G.-J.; deLeeuw, S. W.; Eisenbach, C. D. J. Chem. Phys. 2000, 113, 4. (19) Barett-Jones, D. P.; Leblanc, B.; Hertfort, M.; Moss, T. Science 1994, 264, 1134. (20) Shih, I.-L.; Van, Y.-T. Bioresour. Technol. 2001, 79, 207. (21) Bizzarri, A.; Cametti, C.; Bordi, F. J. Phys. Chem. 1990, 94, 2166. (22) Huh, Y.; Kai, S.; Satoh, M.; Komiyama, J. Macromol. Chem. Phys. 1994, 195, 3813. (23) Wallis, R.; Record, B. Trans. Faraday Soc. 1962, 58, 1251. (24) Holtzer, A.; Hawkins, R. J. Am. Chem. Soc. 1996, 118, 4220. (25) Bianchi, E.; Bicchi, A.; Conio, G.; Ciferri, A. J. Macromol. Sci., A1 1967, 909. (26) Case, D. A. AMBER9.0; University of California: San Francisco, 2006. (27) Jorgensen, W. L.; Chandrasekhar, J.; Madura, J. D.; Impey, R. W.; Klein, M. L. J. Chem. Phys. 1983, 79, 926. (28) Humphrey, W.; Dalke, A.; Schulten, K. J. Mol. Graphics 1996, 14, 33. (29) Fedorov, M. V.; Goodman, J. M.; Schumm, S. J. Am. Chem. Soc. 2009, 131, 10854. (30) Dang, L. X. J. Am. Chem. Soc. 1995, 117, 6954. (31) Doi, M.; Edwards, S. F. The Theory of Polymer Dynamics; Clarendon Press: Oxford, 1986. (32) Finke, J. M.; Jennings, P. A.; Lee, J. C.; Onuchic, J.; Winkler, J. R. Biopolymers 2007, 86, 193. (33) Miettinen, M. S.; Gurtovenko, A. A.; Vattulainen, I.; Karttunen, M. J. Phys. Chem. B 2009, 113, 9226. (34) Berendsen, H. J. C.; Grigera, J. R.; Straatsma, T. P. J. Phys. Chem. 1987, 91, 6269. (35) Smith, D. E.; Dang, L. X. J. Chem. Phys. 1994, 100, 3757. (36) Daub, C. D.; Leung, K.; Luzar, A. J. Phys. Chem. B 2009, 113, 7687. (37) Horinek, D.; Mamatkulov, S.; Netz, R. R. J. Chem. Phys. 2009, 130, 124507. (38) Joung, I. S.; Cheatham, T. E., III J. Phys. Chem. B 2008, 112, 9020. (39) Åqvist, J. J. Phys. Chem. 1990, 94, 8021. (40) Auffinger, P.; Cheatham, T. E., III; Vaiana, A. C. J. Chem. Theory Comput. 2007, 3, 1851. (41) Metzler, R.; Klafter, J. Chem. Phys. Lett. 2000, 321, 238. (42) Barrat, J. L.; Hansen, J.-P. Basic Concepts for Simple and Complex Fluids; University Press: Cambridge, 2003. (43) Fugman, S.; Sokolov, I. M. J. Chem. Phys. 2009, 131, 235104. (44) Alexander-Katz, A.; Wada, H.; Netz, R. R. Phys. ReV. Lett. 2009, 103, 028102. (45) Sangha, A. K.; Keyes, T. J. Phys. Chem. B 2009, 113, 15886. (46) Neusius, T.; Daidone, I.; Sokolov, I. M. Phys. ReV. Lett. 2008, 100, 188103.

JP1010814