Divalent Ion Parameterization Strongly Affects Conformation and

Feb 16, 2016 - Divalent Ion Parameterization Strongly Affects Conformation and Interactions of an Anionic ... Copyright © 2016 American Chemical Soci...
3 downloads 0 Views 4MB Size
Article pubs.acs.org/JPCB

Divalent Ion Parameterization Strongly Affects Conformation and Interactions of an Anionic Biomimetic Polymer Michael D. Daily,†,‡ Marcel D. Baer,*,† and Christopher J. Mundy† †

Physical Science Division, Pacific Northwest National Laboratory, Richland, Washington 99352, United States Department of Biochemistry and Molecular Biology, University of Texas Medical Branch, Galveston, Texas 77555, United States



S Supporting Information *

ABSTRACT: The description of peptides and the use of molecular dynamics simulations to refine structures and investigate the dynamics on an atomistic scale are well developed. Through a consensus in this community over multiple decades, parameters were developed for molecular interactions that only require the sequence of amino-acids and an initial guess for the three-dimensional structure. The recent discovery of peptoids will require a retooling of the currently available interaction potentials in order to have the same level of confidence in the predicted structures and pathways as there is presently in the peptide counterparts. Here we present modeling of peptoids using a combination of ab initio molecular dynamics (AIMD) and atomistic resolution classical force field (FF) to span the relevant time and length scales. To properly account for the dominant forces that stabilize ordered structures of peptoids, namely steric-, electrostatic, and hydrophobic interactions mediated through side chain−side chain interactions in the FF model, those have to be first mapped out using high fidelity atomistic representations. A key feature here is not only to use gas phase quantum chemistry tools, but also account for solvation effects in the condensed phase through AIMD. One major challenge is to elucidate ion binding to charged or polar regions of the peptoid and its concomitant role in the creation of local order. Here, similar to proteins, a specific ion effect is observed suggesting that both the net charge and the precise chemical nature of the ion will need to be described.



INTRODUCTION Background ions in solution are known to impact synthetic outcomes and the hierarchical assembly of complex structures in biology and elsewhere. An example of this is the specific ion effect, which is quantified by the well-known Lyotropic and Hofmeister series and has dramatic consequences for the outcomes of colloidal and protein assembly.1,2 In addition, many proteins require a metal ion to maintain their native structure and/or function.3 While Hofmeister effects are weaker for cations than anions, divalent cations do have weak affinity for the backbone carbonyl oxygen,4,5 but monovalent ions have negligible affinity. Most relevant for this work, there is significant evidence for the importance of mono- and divalent ions such as K+ and Ca2+ in influencing the structures of anionic peptides,6 phosphopeptides,7 peptide hydrogels,8 and allosteric proteins like calmodulin.9 Thus, to develop an accurate theory of the effect of ions on the conformation and structure of macromolecules, it is crucial to understand the local solvation structure of ions, an area which remains theoretically challenging. To this end, a surge of theoretical and computational activity has linked the representation of ion solvation to processes at the air−water interface and colloidal assembly.2,10−18 The simulation-based discovery that halide anions partition to the air−water interface in a reverse Hofmeister series suggested that polarizability was the major microscopic driving force for the partitioning of ions to model hydrophobic interfaces.2 Since © 2016 American Chemical Society

these pioneering simulations, the identification of ion polarizability as the dominant force driving ions to hydrophobic interfaces has been refined.2,10,11,13−15,17−20 Furthermore, ongoing research is converging toward a picture that implicates the subtle balance of hydrophobic cavitation and electrostatics that alter the flexibility of the aqueous first solvation shell.2,14,15,21 Due to this emerging consensus, there is a premium on experimental determination of the intrinsic solvation structure of aqueous ions and the relationship of the local structure to more collective phenomena such as ion-pairing and association.22 Extended X-ray absorption fine structure (EXAFS) measurements have been key to accurately measuring the solvation structure of aqueous ions, resolving the first shell water molecules (position, number, and orientation) as well as effects of concentration due to ion pairing in electrolytes. Almost universally, quantum mechanical forms of molecular interaction based in density functional theory (DFT) have outperformed their classical counterparts in reproducing experimentally observed EXAFS spectra.23−31 Advances in the efficiency of DFT methods have allowed researchers to study more complex processes beyond single ion solvation such as ion pairing.5,30,32−36 These improved DFT Received: December 15, 2015 Revised: February 5, 2016 Published: February 16, 2016 2198

DOI: 10.1021/acs.jpcb.5b12277 J. Phys. Chem. B 2016, 120, 2198−2208

Article

The Journal of Physical Chemistry B

Most importantly, we find that Ca2+ - model peptoid PMFs differ significantly between DFT and classical simulations. We import these findings from DFT and modify the classical FF to improve agreement with the DFT-based PMFs using an approach similar to the ECC method discussed above. We find that the conformational ensemble of the anionic biomimetic polymer differs dramatically depending on the form of the ion-acetate and ion-carbonyl interactions. Our studies suggest that improved free energy estimates for monoand divalent ion interactions with backbone carbonyls and negatively charged side-chains could lead to improved classical models for systems whose structures depend on ion-mediated interactions.

methods have produced a picture of these fundamental process that is distinct from the description obtained from classical empirical potentials. For example, a full ab initio treatment of pairing between cations Na+ and Ca2+ and peptide backbone analogue n-methylacetamide (NMA) explained the experimentally observed infrared (IR) shift of the amide I band as a function of electrolyte.5 The DFT-produced potentials of mean force (PMF) and local aqueous structure for ion binding to NMA were key to understanding how the specific ion effect is manifest in IR spectra and thus on structure of model proteins. In particular, classical simulations with Ca2+ often perform poorly by comparison to quantum mechanics based simulations because this ion strongly polarizes neighboring water molecules.37 Recently, Leontyev et al. have shown that classical molecular dynamics (MD) can approximately account for polarization effects using an “electronic continuum correction” (ECC) that scales ion charges by (1/ϵel) , or about 0.75, where ϵel is the electronic part of the solvent dielectric constant.38 Jungwirth and co-workers showed that classical simulations of CaCl2 using a revised version of ECC predict significantly better pair distribution functions than were obtained with standard empirical interaction potentials.39,40 In addition, Kohagen et al. showed that classical MD simulations of calmodulin with a revised ECC reproduced experimental data much better (calcium binding affinity and the structure of the Ca2+ binding sites) than standard force fields (FF).41 While the utility of ECC as a general approach to incorporate manybody effects remains to be seen, there is strong evidence for its utility for systematic testing and correcting of ion−ion interactions to yield the experimentally determined collective response.22 In the realm of bioinspired materials, there remains interest per se in engineering rational foldamer designs,42−44 but this perhaps has been surpassed due to the recent interest in bioinspired self-assembled materials.45,46 For example, peptoids, which are N-substituted glycines, are a highly useful biomimetic polymer (BMP) system due to their high thermal and chemical stability,47,48 resistance to degradation, relative lack of toxicity, and ease of synthesis.49−52 In addition, while synthetic polymers have not so far achieved the same level of structural and functional complexity and specificity as proteins, BMPs have a greater capacity to form extended periodic nanostructures.45 Since peptoids bear a side chain on the backbone nitrogen rather than on the Cα, they lack backbone chirality and the ability to form backbone−backbone hydrogen bonds like peptides.53,54 Our research presented in here shows, that the non-hydrogen-bonded carbonyl groups give peptoids an enhanced divalent ion sensitivity by comparison to proteins. Thus, for this system, it will be imperative to have a deep understanding of how the solution conditions impact the hierarchical assembly of materials such as membranes, fibers, and pores for large-scale molecular separation or catalysis. In this work, we aim to better understand ion-mediated assembly and folding in these peptoids and the competition with hydrophobic effects. We systematically investigate how the conformational populations of an anionic BMP with a welldefined sequence depend on the precise treatment of the background electrolyte comprised of CaCl2 or KCl. Detailed potential of mean force (PMF) calculations of Ca2+ and K+ ions with the model peptoid dimethylacetamide (DMA) and the acetate anion inform our understanding of the nascent interaction.



METHODS Parameters for Peptoids and Analogues. Since peptoid topologies and parameters are not present in standard MD FF, we generated them using the generalized Amber force field (GAFF) and the Antechamber program.55,56 The AM1-BCC model57,58 was used to generate partial charges as used in previous studies in implicit solvent.59 GAFF-based models have been shown to correctly predict peptoid crystal structures.60 Recent studies that focus to improve the backbone potentials based on the AMBER61 and CHARMM force field46,62 highlight the differences in conformational space compared to peptides. For this study we did not reparameterize the backbone potentials but focused on the interaction for ionic side chains with electrolytes. GAFF parameters were converted into GROMACS-compatible topologies using the AcPype progam, 63 and simulations were carried out using AMBER0364 with the added peptoid parameters. Since AMBER NaCl and KCl are prone to aggregating into salt crystals around biomolecules,65 we use the Dang-Smith potentials for potassium ions instead.66 Dang et al. did not modify the chloride parameters. To scale the charge of carboxylate groups by 0.75, the charges on the three carboxylate atoms (CD, OE1 and OE2) and the adjacent methylene group (CG, HG1 and HG2) in the case of Nce) are uniformly scaled to add 0.25 to the total charge of the residue so that the total charge is −0.75, see Table 1 and SI for details. Table 1. Original and Scaled Charges on the Nce Residue atom

CG

HG1/2

CD

OE1/2

charge unscaled charge scaled

−0.229400 −0.171762

0.032700 0.024484

0.858602 0.642873

−0.844801 −0.632539

We include the topologies and FF parameters for all peptoids and other molecules in this work as an attachment in the Supporting Information. First, peptoids were constructed in extended conformations with all-trans backbone amides using an in-house python script. Peptoids are capped on the Nterminus by an acetyl group and on the C-terminus by a N,Ndimethyl group. These were then equilibrated in implicit solvent as described below. For peptoid monomer simulations, we solvated the starting structure in a 15 Å thick layer of SPC water67 in a dodecahedral box and then added the desired number/concentration of ion(s). To simulate peptoid dimers, the corresponding monomer simulations were clustered to identify favored conformations as described below, and then two copies of the top or second-ranked cluster were placed 20 Å apart, and the box was expanded by an additional 10 Å in 2199

DOI: 10.1021/acs.jpcb.5b12277 J. Phys. Chem. B 2016, 120, 2198−2208

Article

The Journal of Physical Chemistry B

to 5.0 Å equally spaced by 0.2 Å for the O···cation and C··· cation, respectively. A harmonic umbrella potentials of the form V(r) = k(r0 − r)2 with a force constant k of 9.56 Å−2 kcal mol−1 was used. In each umbrella window a trajectory of at least 40 ps was collected after 5 ps of equilibration. The weighted histogram analysis method (WHAM) was employed to extract a free energy profile from these histograms.74 For calculating potentials of mean force between carboxylate groups in the presence of ions, we used the peptoid analogue MCX (see Figure 1). Two MCX molecules were spaced 10.0 Å apart, and the system was solvated with a 15 Å water layer. The desired concentration of ions was added and the simulation was carried out as described above for 2 μs. To calculate the PMF, the distribution of the center-of-mass distance between the two carboxylate groups (CG, OD1, OD2) was calculated with 0.5 Å bins. The densities were averaged over 100 ns blocks, and the free energies were calculated by ΔGi = −kBT(ln(P(ri)) − 2 ln(4π ri)), where ri is the midpoint of bin i, P(ri) is the probability density for bin i, 2 ln(4π ri) the volume correction factor, and ΔGi is the ΔG at ri. Clustering. To identify favored conformations from the simulation, we performed k-means clustering75 of the O (1000) ensemble structures in R,76 using the root-mean-square distance of the backbone atoms (N, CA, C, and O) as a distance metric. All structures in the ensemble (after equilibration time) are used for clustering. To optimize the number of clusters k, we clustered using a range of k and identified the “knee” in a plot of F(k) vs k, where F(k) is a “goodness of clustering” quality metric. Ten random starting assignments were used at each k. First, we quantify the quality of each individual cluster using the “internal RMS”:

each dimension before adding water and/or ions. The use of prefolded monomers can be thought of as reflecting the dilute limit in which monomers are likely to become structured before binding to each other. For implicit solvent simulations, we used the GBSA implicit solvent model with the Onufriev/Bashford/Case algorithm68 for calculating Born radii and a dielectric of 78.3, and infinite van der Waals and Coulomb cutoffs. We used nonperiodic NVT simulations with velocity-rescale temperature coupling and a 0.2 ps coupling time.69 Simulation Protocol for Classical Simulations. MD simulations were performed with GROMACS 4.6.4.70 We used 10 Å cutoffs for van der Waals and Coulombic interactions and the particle-mesh Ewald method71 for long-range electrostatic interactions. The simulations were performed in the NPT ensemble with a Parrinello−Rahman barostat72 with a 1 ps coupling time at 1 bar. A Nosé−Hoover thermostat73 at 300 K with 0.2 ps coupling time was applied separately to the polymer and to the solution (including ions). The simulations were performed for variable times with a time step of 2 fs. We perform simulations with two parameter sets that differ in both the charge on the Ca2+ and the partial charges on the carboxylate groups in the BMP. We will refer to the standard FF with a superscript s, and the revised parameters with a superscript r. For analyses, we discard the first 500 ns of Cas Cls2 simulations and the first 100 ns of other simulations as equilibration time. We use the larger cutoff for Cas Cl2s simulations because the number of ion bridges typically takes about this long to equilibrate (not shown). The number equilibrates much more quickly in Car Clr2 due to the weaker Car-carboxylate interactions. The N-methyl,N-methylcarboxyl-acetamide (MCX, see Figure 1a) simulation cell consisted of one solute molecule in 206 water molecules and a single Ca2+, Na+, and K+ ion in a cubic box with the edge length of 20.00 Å. The free energy profiles for the C···cation and CO···cation interaction were obtained by umbrella sampling in two independent simulations. Sampling windows for distances range from 2.0 to 6.0 Å and 2.8

Ns

⟨∑ ∑ dij 2⟩

RMSint =

i=1 j>i

(1)

where RMSint is the internal RMS, Ns is the number of structures in the cluster, and dij is the root-mean-square backbone distance between structures i and j. Then, for F(k), we calculate the “order ratio” as follows: k

F (k ) =

∑i = 1 NC i i k

∑i = 1 Ni

(2)

where Ni is the number of structures in cluster i and Ci is a Boolean that is 1 if RMSint for cluster i is less than some arbitrary cutoff and 0 otherwise. For the clustering calculations in this work, we use an RMSint cutoff of 2.5 Å unless otherwise specified. Ab Initio Simulations. Born−Oppenheimer ab initio molecular dynamics simulations within the NVT (at 300 K) ensemble using periodic boundary conditions are perfomed within the CP2K simulation suite (http:www.cp2k.org) containing the QuickStep module for the DFT calculations.77 We followed a similar protocol as in ref 5, using a double-ζ basis set that has been optimized for the condensed phase78 in conjunction with GTH pseudopotentials79 using a 400 Ry cutoff for the auxiliary plane wave basis. A Nosé−Hoover thermostat was attached to every degree of freedom to ensure equilibration.80 The Becke exchange81 and correlation due to Lee, Yang, and Parr (LYP)82 is utilized in addition to the dispersion correction (D2) put forth by Grimme83 with a 40 Å, cutoff. The DMA simulation cell consisted of one solute molecule in 107 water molecules and a single Ca2+, Na+ ion in a

Figure 1. Peptoid analogue used in simulating ion-mediated interactions. (a) force field compound representing charged side chain and backbone carbonyl, (b) acetate, and (c) N,N-dimethylacetamide (DMA) used to study cation peptoid interaction using DFT (see Methods for details). 2200

DOI: 10.1021/acs.jpcb.5b12277 J. Phys. Chem. B 2016, 120, 2198−2208

Article

The Journal of Physical Chemistry B

referring to calcium and chloride with set 3 (revised) parameters as Car and Clr, respectively.

cubic box with the edge length of 14.96 Å. The free energy profiles for the CO···cation interaction were obtained by umbrella sampling. Sampling windows for the O···cation distance (r0) ranging from 2.2 to 5.6 Å were equally spaced by 0.2 Å employing harmonic umbrella potentials of the from V(r) = k(r0 − r)2 with a force constant k of 9.56 Å−2 kcal mol−1. To ensure sufficient sampling in the barrier region for the Ca2+ additional windows with stiffer force constants were added ranging from 2.8 to 4.1 Å equally spaced 0.1 Å apart with a force constant of 19.12 Å−2 kcal mol−1. In each umbrella window a trajectory of at least 25 ps was collected after 5 ps of equilibration. The acetate simulation cell consisted of one solute molecule in 95 water molecules and a single Ca2+, Na+, and K+ ion in a cubic box with the edge length of 14.40 Å. The free energy profiles for the C···cation interaction were obtained by umbrella sampling. Sampling windows for the O···cation distance ranging from 2.4 to 5.6 Å were equally spaced by 0.2 Å employing harmonic umbrella potentials of the form V(r) = k(r0 − r)2 with a force constant k of 9.56 Å kcal mol−1. To ensure sufficient sampling in the barrier region for the Ca2+ and Na+ additional windows with stiffer force constants were added ranging from 3.6 to 4.6 Å equally spaced 0.1 Å apart with a force constant of 33.46 Å−2 kcal mol−1. In each umbrella window a trajectory of at least 40 ps was collected after 5 ps of equilibration. The weighted histogram analysis method (WHAM) was employed to extract a free energy profile from these histograms.74 Backbone and Side-Chain Interactions with Ions. Here we study the ion induced structuring of peptoids with alternating hydrophobic and anionic side chains. As has been shown for proteins, the specific effects of cations on the structure are primarily due to interactions with anionic side chains and backbone carbonyls. The basic building unit of the anionic residue that is also studied experimentally is shown in Figure 1. We start by comparing binding of mono and divalent cations with the model residue. Since peptoid topologies and parameters are not present in standard MD FFs, we developed a new set (see Methods and SI for details). For the DFT simulations we investigate the binding to the neutral and anionic part of the residue using two independent simulations. To study the binding of mono- and divalent cations to the peptoid backbone carbonyl we used N,N-dimethylacetamide (DMA) as model system similar to the NMA in a recent study for the peptide backbone.5 As a model for the anionic side chain we use the acetate anion as suggested in ref 35. The largest ion induced effect on the conformation and assembly for the charged model peptoids is expected for the ion interaction with the anionic side chain, as shown in the recent literature capturing both the local solvation structure and longrange collective response within a single parametrization for Ca2+.22,39,84 As there is no experimental data available that we can use to compare the strength of ion binding to the peptoids we use DFT as reference. For three different sets of parameters we compare structural motifs that are formed for the model peptoids and their residence times. Set 1 is AMBER03 that is used in conjuction with GAFF. Following the approach of refs 39−41, we use a charge scaling approach reducing the charge on the Ca2+ by 25%, but keeping the Lennard-Jones parameters unchanged, denoted set 2. Set 3 a charge scaling is also performed on anionic solutes, see SI for details. For reference, we will be referring to calcium and chloride with standard AMBER03 parameters (set 1) as Cas and Cls, respectively, for the rest of the manuscript. We will also be



RESULTS AND DISCUSSION We start our discussion with comparing ion binding for the acetate with K+ and Ca2+. In Figure 2a the potential of mean

Figure 2. PMF of calcium (panel a), sodium (panel b), or potassium (panel c) and carboxylate as a function of the distance between the calcium cation and the carbon atom of the carboxylate group of the acetate and N-acetyl-N-methylglycinate for BLYP-D and MM, respectively. The curves have been shifted such that their minimum (either bidentate or monodentate) is aligned at 0 kcal/mol. Arrows up and down in panel a indicate the position of the bidentate and monodentate contact ion pair, respectively.

force (PMF) for the binding of Ca2+ to the acetate is shown using the ion···C (carbon of the carboxylate group) distance as reaction coordinate following a similar approach as in ref 35 The PMFs in ref 35 and our current results are in qualitative agreement given the large uncertainties in the earlier study and the different functional. Here, we find that both mono- and bidentate motifs (see Figure 3) are equally stable for the contact ion-pair (CIP), which is about 1 kcal/mol more stable than the solvent separated ion-pair (SSIP). The occurrence of both binding motifs is known to play a role in the structure and function of metalloproteins.85 The barrier for dissociation for our current simulations and protocol are in good agreement 2201

DOI: 10.1021/acs.jpcb.5b12277 J. Phys. Chem. B 2016, 120, 2198−2208

Article

The Journal of Physical Chemistry B

Figure 3. Representative structures along DFT simulations for the acetate−Ca2+ PMF corresponding to (a) bidentate contact ion pair, (b) monodentate contact ion pair, and (c) solvent separated ion pair. Thin gray bond indicate oxygen atom that are within the first solvation shell (distance smaller than 2.6 Å) of the Ca2+-ion. For these structures corresponding to the minima on the free-energy surface shown in Figure 2a the calcium ion is predominantly found with a 7-fold coordination similar to the bulk solvation structure as shown in ref 87.

with a previous study, namely 2.6 and 3.2 kcal/mol, respectively.35 The most pronounced difference is the stability of the SSIP, which is almost 1.8 kcal/mol more stable in our current simulation. For structures corresponding to the minima of the free-energy profile shown in Figure 2a the calcium ion is predominantly found with a 7-fold coordination similar to the bulk solvation structure as shown in refs 29, 86, and 87. Additionally we find structures with both six and eight oxygens in the first solvation shell are sampled. These states are accessible in coordination number free-energy landscape for Ca2+ at room temperature.87,88 Comparing the DFT and AMBER03 simulation (set 1) we find qualitative differences for the stability of the SSIP that is found at distances of 4.5 and 5.0 Å, respectively (Figure 2a). This results in a large population shift toward the contact ion pair for set 1 with an increase for the barrier of dissociation by a factor of 2 over the DFT of 5.8 kcal/mol. Scaling the charge on the cation (set 2) leads to much better agreement between the DFT and classical FF PMFs. The relative stability between the SSIP and CIP is equal within the uncertainty and the barrier between the CIP and SSIP moiety is also roughly equal. There is, however, a noticeable shift of the PMF obtained from set 2 to larger distances. Scaling the charge on the solute (set 3) corrects for this structural change. For the case of the two monovalent cations (potassium and sodium) considered herein, good agreement is already observed in the PMFs between the DFT and the standard classical FFs as shown in Figure 2b,c. For K+ the PMFs for the DFT and the classical FF are the same within the errorbars. The Na+ shows the bidentate CIP 0.5 kcal/mol more stable for the classical FF than the global minimum mono dentate configuration for the DFT. Nevertheless, the barrier from the SSIP to CIP is well reproduced between classical and DFT based interaction potentials. For both K+ and Na+ the CIP is more stable by 1 kcal/mol than the SSIP. The binding to the backbone carbonyl of peptoids is studied using DMA as the model molecular fragment. The PMF of ion binding is shown in Figure 4a,b for Ca2+ and Na+. Not surprisingly, the predicted binding to the carbonyl is similar to that observed in a previous NMA study,5 the peptide backbone analogue. The CIP between Ca2+ and the peptoid carbonyl is 1.5 kcal/mol more stable than the SSIP and differs from the NMA result by roughly 1 kcal/mol.5 The same trend is seen for the Na+ binding, where the CIP and SSIP have the same freeenergy. Interestingly, the barrier for association is 3 kcal/mol for both Na+ and K+ ions. Comparing the PMF for ion-binding for the two sets of Ca2+ parameters, it can be seen, that the

Figure 4. PMF for pairing between DMA and Ca2+ (panel a). Comparing BLYP-D (black) to the original AMBER parameter set (red) and the scaled charge approach (green). Panel b shows freeenergy of DMA and Na+ (BLYP-D in black and AMBER in blue) and K+ (Dang-Smith (DS) in green) paring. See Methods section for details.

original parametrization (set 1) gives almost an identical PMF to the DFT and where only the minimum of the CIP is shifted by about 0.2 Å. Scaling the charge leads to a reduction in the barrier between the CIP and SSIP states from 5 to 3 kcal/mol, but keeps the relative stability between the two states. This may have an effect on the observed kinetics of ion binding, however, herein we are primarily interested in getting the relative populations of the CIP and SSIP for both the carbonyl and acetate groups with respect to Ca2+. Once again, we find the unmodified classical FF for the monovalent ion−carbonyl PMF is in good agreement with that obtained using DFT. Thus, we only need to modify the Ca2+ to best reproduce populations of both ion−carbonyl and ion−acetate populations that are consistent with DFT. Effect of Model Choice on Apparent Free Energy of Ion-Mediated Carboxylate Interactions. To directly compute the free energetics of cation-mediated interactions between carboxylates, we first calculate the PMFs between two methyl-carboxylate side chain (Figure 5) as a function of both the electrolyte and ionic strength as described in the Methods. Figure 5a shows that with 0.2 M KCl, the effective interactions between carboxylates are purely repulsive. By contrast, at Cas concentrations up to 0.10 M, the carboxylates interact at a distance of about 5.0 Å with an effective attractive free energy of 1−1.5 kcal/mol. This is consistent with a direct −CO2−···Cas ··· −2OC- bridge. The 0.03 M Cas simulation uses only one Cas, just enough to neutralize the two carboxylates. Surprisingly, at 0.50 M Cas Cls2, the effective free-energy of association between −CO2− is less than 0.5 kcal/mol. This occurs because strong interaction between Cas, essentially saturates all available −CO2− at high concentration and there probability of finding 2202

DOI: 10.1021/acs.jpcb.5b12277 J. Phys. Chem. B 2016, 120, 2198−2208

Article

The Journal of Physical Chemistry B

Figure 5. Effective potentials between carboxylates in the presence of various ions. The PMFs are calculated from unbiased atomistic simulations of 2 MCX as described in the Methods. The errors are standard error of the mean from averaging over 100 ns blocks. PMFs with 0.20 M KCl (black) and (a) varying concentrations of Cas PMFs and (b) with varying concentrations of Car.

a −CO2−···Cas ··· −2OC- bridge is low. Figure 5b shows that effective inter--CO2− interactions are much weaker in the presence of Car than with Cas. At 0.03 M Car, the interaction is nearly zero, while at 0.50 M it is only 0.39 kcal/mol at 5.8 Å. However, the seemingly unphysical inverse concentration dependence seen with Cas is not observed with Car. These results suggest that the scaling of the Ca2+ charge reduces a strong ion-mediated intercarboxylate interaction to a freely reversible one. Moreover, there are no significant ion bridges between the carbonyls of two MCX observed for either FF parameters for Ca2+. Effect of Ions and Model Choice on Peptoid Backbone Conformations. We first examine the contribution to folding from cation-carbonyl interactions using sar12, a 12-mer of sarcosine (Figure 6a). Figure 7 shows that in pure water as well as in 1.0 M KCl, sar12 Rgyr averages about 8.5 Å. By contrast, 0.10 M Cas Cls2 slightly shifts the conformational distribution toward more compact conformations, and 0.50 M Cas Cls2 greatly increases

Figure 7. Representative structure for the calcium-coil (panel a) and concentration dependence of the radius of gyration (Rgyr) for a neutral model peptoid (all-trans sar12) (panel b). Radius of Rgyr is calculated from all Cα atoms of all residues. Black, green, blue with open symbols, and red with open symbols indicate no salt, 1.0 M KCl, 0.10 M Cas Cls2 and 0.50 M Car Clr2, respectively. Comparison of Rgyr distribution using Cas (open symbols) vs Car (closed symbols). Salt-free and 1.0 M KCl simulations were run for 0.5−1 μs, and CaCl2 simulations were run for 1−2 μs. The errors, where given, are standard error of the mean from averaging over 50 ns blocks.

the population of compact conformations (Rgyr < 7.0 Å) by comparison to water and 1.0 M KCl. Figure 7b illustrates a “calcium coil” structure of sar12 with Rgyr = 5.6 Å that rarely occurs in 0.5 M Cas Cls2. The backbone coils around two central calcium ions, and the four carbonyl oxygens are within 4.0 Å of

Figure 6. Model 12-mer peptoids with identical anionic side chains and increasingly hydrophobic apolar side chains. (a) sar12, (b) (Nce-Net)6, (c) (Nce-Nib)6, and (d) (Nce-Ncp)6. 2203

DOI: 10.1021/acs.jpcb.5b12277 J. Phys. Chem. B 2016, 120, 2198−2208

Article

The Journal of Physical Chemistry B Table 2. Ion-Binding Statistics for Three Peptoids of Varying Hydrophobicitya average no. bound to Ca system s

(Nce-Net)6 Ca (Nce-Net)6 Car (Nce-Nib)6 Cas (Nce-Nib)6 Car (Nce-Ncp)6 Cas (Nce-Ncp)6 Car

average no. Ca bound to

carbonyl

-CO2−

carbonyl

-CO2−

total

no. ion bridges

favored structure

1.35 7.58 1.21 7.20 1.65 6.11

5.96 1.87 5.95 1.81 5.84 1.17

0.64 2.01 1.07 1.89 0.93 1.30

5.14 2.19 4.95 1.81 5.21 1.18

5.26 2.94 4.99 2.90 5.69 2.32

7.66 0.23 7.86 0.78 4.37 0.04

disordered with many ion bridges calcium coil disordered bb with many ion bridges ion bridges mix of Ca coil and less ordered structures ion bridges calcium coil

Nib indicates N-isobutyl glycine, and Ncp indicates N-ethylphenylchloryl glycine. The first 500 ns of each simulation are excluded from statistical calculations. a

either Cas. Based on Coulombic forces alone it follows that Car will interact more weakly with the peptoid backbone than Cas, leading to a reduced compacting effect. However, Figure 7a shows that with 0.1 M Car Clr2, the conformational distribution is slightly more compact than with 0.1 M Cas Cl. At 0.5 M, the additional compacting effect of Car and Cas is even larger, and the dominant peak is the calcium coil at 6.0 Å. Figure 7b illustrates the cluster with highest population from the 0.5 M Car Clr2 simulation, which is a more compact calcium coil than that observed with Cas Cls2. Nine carbonyl oxygens are within 4.0 Å of either Car. Effect on Structure of Anionic Polymers with Different Hydrophobic Side Chains. To connect the energetics of ion-mediated interactions to peptoid structure, we examine how the conformation of anionic, amphiphilic peptoid monomers respond to calcium ions with different parameters. We also examine how the ion-mediated interactions cooperate or compete with hydrophobic interactions to govern structuring. These model peptoids have sequences (NceNhp)6, where Nce is an N-ethylcarboxylate residue and Nhp is a side chain of variable hydrophobicity. For the different peptoids (shown in Figure 6), Table 2 shows how peptoid/ Ca2+ binding and ion-mediated interactions vary among models. The first 500 ns of the simulations are excluded from these statistical calculations. First, to examine ion-mediated assembly in isolation from hydrophobic effects, we examine (Nce-Net)6, which contains weakly hydrophobic N-ethyl side chains. For this peptoid in 0.50 M Cas Cls2 on average, 1.4/13 carbonyls bind to less than one Cas, all six −CO2− bind to about 5.1 Cas, and there are about 7.7 ion bridges. By contrast, in 0.50 M Car Clr2, an average of 7.6/13 carbonyls bind to about 2 Car, while only 1.87/6 −CO2− bind to about 2.2 Car, and the ion bridges are nearly absent. The carbonyl statistics for this system are consistent with the formation of a calcium coil like that illustrated in Figure 7b. In addition, the ion bridge results are consistent with the much weaker ion bridges observed for model compounds in Car Clr2 than in Cas Cls2 (Figure 5). These results show that (Nce-Net)6 switches from a structure dominated by Cas/-CO2− interactions and ion bridges in Cas Cls2 to a structure dominated by Car/carbonyl interactions in Car Clr2. Cluster analysis shows that in Cas Cls2, with an internal rmsd cutoff of 2.5 Å, there are two ordered clusters comprising 98% of the ensemble after 500 ns. Both clusters have the same ion-binding profile and similar structures; Figure 8a shows that cluster 2 has many Cascarboxylate interactions and ion bridges. By contrast, Figure 8b shows that for in Car Clr2, 99% of the ensemble after 100 ns forms a calcium coil similar to the sar12 structure in Figure 7b.

Figure 8. Variation of different anionic peptoids’ conformations with calcium model. (a) Cluster 2 from simulation of (Nce-Net)6 with 0.50 M Cas Cls2, with bound ions included. (b) Cluster 1 from simulation of (Nce-Net)6 with 0.50 M Car Clr2. (c) Cluster 3 from simulation of (Nce-Nib)6 with 0.50 M Car Clr2. (d) Cluster 2 from simulation of (Nce-Ncp)6 with 0.50 M Car Clr2.

Next, we examine (Nce-Nib)6, for which the hydrophobic isobutyl side chains are of similar size to the anionic Nce side chains. Considering this, it is plausible that a stable hydrophobic core can form. Table 1 shows that for this peptoid, the ion-binding properties for each Ca model are very similar to (Nce-Net)6, except that slightly more ion bridges form using Car. The carbonyl-Car binding statistics are consistent with the Car Clr2 ensemble being dominated by the calcium coil motif. In Cas Cls2, (Nce-Nib)6 forms similar structures to those observed for (Nce-Net)6 in Cas Cls2. In Car Clr2, 28% of the ensemble after 100 ns is disordered, while the top two ordered clusters comprise 46% of the ensemble and adopt calcium coil structures. The remaining two ordered clusters (comprising 26%) form less ordered calcium coils, as illustrated in Figure 8c for cluster 3. Finally, the data unexpectedly show that most Nib side chains remain solvent-exposed using either Cas Cls2 or Car Clr2. Finally, we examine (Nce-Ncp)6, where the Ncp side chain is highly hydrophobic and twice the size of Nib and Nce. Thus, even if this peptoid were to form a hydrophobic core, much of the hydrophobic surface area would still be solvent-exposed. Table 2 shows that for this peptoid, the ion-binding properties 2204

DOI: 10.1021/acs.jpcb.5b12277 J. Phys. Chem. B 2016, 120, 2198−2208

Article

The Journal of Physical Chemistry B for each calcium model are very similar to the other peptoids, except that this peptoid forms fewer ion bridges in Cas Cls2 than either of the other peptoids. (Nce-Ncp)6 forms similar structures in Cas Cls2 as (Nce-Net)6, except that some stacking of Ncp side chains occurs. In Car Clr2, 82% of the ensemble (after 100 ns) is in one of two ordered clusters that form a calcium coil with some Ncp stacking, as illustrated in Figure 8d. In summary, these results show that in anionic peptoids, preferred backbone conformation differs greatly between Cas Cls2 and Car Clr2, but the preferred structures vary little with the hydrophobicity of the Nhp side chain. Considering this and that the (Nce-Nib)6 and (Nce-Ncp)6 simulations show substantial exposed hydrophobic surface area, our data suggest that peptoid/cation interactions are more important than hydrophobic interactions to the folding of monomers. Effect on Dimerization of Anionic Polymers with Different Hydrophobic Side Chains. To illuminate initial stages of ion-mediated peptoid assembly, we examine the calcium model-dependence of (Nce-Net)6 dimerization. To mimic the dilute limit, which corresponds most closely to reality, we start these simulations from “pre-folded” top clusters of the corresponding monomer simulations. For the (Nce-Net)6 dimer in 0.50 M Cas Cls2, on average, 5.6/26 carbonyls bind to about 2.3 Cas, see Table 3. On

Figure 9. Dependence of model anionic peptoid (Nce-Net)6 dimer on concentrations of different ions. (a) Cluster 1 from the (Nce-Net)6 dimer simulation in 0.5 M Cas Cls2; this simulation was started from two folded monomers (cluster 1 of the monomer Car simulation). (b) Cluster 1 from (Nce-Net)6 dimer simulation in 0.5 M Cas Cls2, starting from two folded monomers (cluster 1 of the monomer Ca r simulation). (c) Cluster 1 of the (Nce-Ncp)6 dimer simulation in 0.5 M Cas Cls2. (d) Cluster 1 of the (Nce-Ncp)6 dimer simulation in 0.5 M Car Clr2. The ordered cluster cutoff was increased to 5.0 Å for (Nce-Net)6 dimers.

Table 3. Ion-Binding Statistics for Dimers of Two Peptoids of Varying Hydrophobicity average no. bound to Ca system (NceNet)6,Cas (NceNet)6,Car (NceNcp)6,Cas (NceNcp)6,Car

with the monomer. These are similar to the numbers for (NceNet)6. Figure 9c,d show that using either Cas or Car, a small hydrophobic core forms at the dimer interface. With Cas, the peptoid occupies one ordered cluster comprising 89% of the ensemble after 500 ns. Panel C shows the cations mostly bind to the exterior of the dimer where the carboxylates are located, while the Ncp residues cluster toward the center. With Car, a single ordered cluster has a lifetime of about800 ns and occupies 45% of the ensemble after 100 ns. This differs strikingly from (Nce-Net)6 in 0.50 M Car Clr2, where the monomers are usually separate. Figure 9d shows that by contrast to the structure in Cas Cls2, in the Car Clr2 structure, some cations and carboxylates are present at the dimer interface. That is, the Car ions stabilize the −CO2− ions and allow them to participate in an interface stabilized by both hydrophobic and ion-mediated interactions. We quantify the interaction between monomers by Boltzmann-inverting the distributions of center-of-mass distance between monomers (rCM) into PMFs. For reference, the neutral sar12 peptide has no free energy minimum corresponding to a dimer. With Cas, the (Nce-Net)6 monomers interact with a free energy of about 1.0 ± 0.1 kcal/mol in a broad but rugged free energy well ranging from 13 to 22 Å. With Car, the monomers interact reversibly with a smaller free energy of −0.67 ± 0.16 kca/mol at rCM 14−16 Å. It is conceivable that higher apparent free energies of interaction might be observed if prefolded monomers were not used, but sampling the vast combined folding/binding space is beyond the scope of this work. The attractive energy between two monomers is about 50% more than that between two carboxylates in the presence of 0.5 M Car Clr2 (Figure 5a). (Nce-Ncp)6 has a dimerization free energy of −3.4 ± 0.14 kcal/mol in 0.50 M Car Clr2.

average no. Ca bound to

carbonyl

-CO2−

carbonyl

-CO2−

total

no. ion bridges

5.7

11.9

2.3

9.7

10.0

17.5

14.5

3.8

3.8

4.2

5.7

0.6

3.0

11.9

2.9

10.9

11.0

16.8

15.7

3.8

3.8

3.8

5.4

1.6

average, nearly all 12 −CO2− bind to about 9.7 Cas, and there are about 17.5 ion bridges. By contrast, in 0.50 M Car Clr2, an average of 14.5/26 carbonyls bind to about 3.8 Car, while only 3.8/12 −CO−2 bind to about 4.2 Car, and there are very few ion bridges, just like with the (Nce-Net)6 monomer. For cluster analysis of these peptoid dimers, we increased the RMSint cutoff for an ordered cluster to 5.0 Å. In Cas Cls2, 45% of the (NceNet)6 simulation ensemble after 500 ns occupies one of two ordered clusters that have many Cas-carboxylate interactions and ion bridges; Figure 9.a illustrates cluster 1. In Car Clr2, the centers of mass of the monomers are only within 20 Å for 12.5% of the simulation ensemble. Thus, we only cluster this subensemble. About 75% of the subensemble occupies one of 10 ordered clusters, all of which maintain the coil structures of the monomers. Figure 9.b illustrates the top-ranked cluster. In the other 9 clusters, the relative orientations of the two monomers are essentially randomly distributed. Table 3 shows that for the (Nce-Ncp)6 dimer in 0.50 M Cas Cls2, on average, 3.0/26 carbonyls bind to about 2.9 Cas. On average, nearly all 12 −CO2− bind to about 10.9 Cas, and there are about 16.8 ion bridges. By contrast, in 0.50 M Car Clr2, an average of 15.7/26 carbonyls bind to about 3.8 Car, while only 3.8/12 −CO2− bind to about 3.8 Car, and there are very few ion bridges, just like



CONCLUSION In summary, we have developed a model that improves the understanding of how hydrophobic effects and ion-mediated interactions cooperate to drive assembly and folding in 2205

DOI: 10.1021/acs.jpcb.5b12277 J. Phys. Chem. B 2016, 120, 2198−2208

Article

The Journal of Physical Chemistry B

national laboratory operated by Battelle for the U.S. Department of Energy. M.D.B. acknowledges support from U.S. Department of Energy, Office of Science, Office of Basic Energy Sciences, Division of Material Sciences & Engineering. C.J.M. acknowledges support from U.S. Department of Energy, Office of Science, Office of Basic Energy Sciences, Division of Chemical Sciences, Geosciences & Biosciences.

peptoids that are used to construct bioinspired materials. Our major finding is the dramatic variation of both internal and collective structure of the BMP based on the precise form of the ion mediated interaction between Ca2+ and the BMP. Specifically, DFT-based PMFs for ion-binding with model peptoid fragments (dimethylacetamide (DMA) and the acetate anion) drive the development of a modified classical FF in which we scale the partial charges of the ions using ideas that are consistent with the ECC method.38 Our findings suggest that by scaling the partial charges for the divalent calcium ion and negatively charged peptoid side chains, we can devise a revised classical FF that much better reproduces the DFT-based PMFs. The work also demonstrates the usefulness of the scaled-charge model for proteins and other polymers and underscores the broad importance of accounting for polarization around divalent ions in solutions of polymers. We find significant differences in the conformational states between the standard and scaled-charge models. Specifically, with the standard model, Ca2+-mediated interactions between the side-chain carboxylates dominates the conformational landscape. However, with the scaled-charge model, the peptoid prefers to form structures with many Car- backbone carbonyl interactions. Such “calcium coil” conformations in the revised FF are favored regardless of the hydrophobicity of the hydrophobic side-chain. Regarding dimers, using the Cas model, monomers strongly interact to make stable structures via ion-bridges that interconvert slowly. By contrast, with Car, the dimers are freely reversible for (Nce-Net)6. For (Nce-Ncp)6 stable dimers form using either Cas or Car, though the latter is more believable given the other evidence presented herein that Car performs better than Cas. The (Nce-Ncp)6 dimer in Car Clr2 is stabilized by both hydrophobic and ion-mediated interactions (Figure 9d), which presents a preliminary picture of how the early stages of peptoid assembly proceed in solution. The disparity of predicted structures and their stability based on the different underlying models could allow for experimental validation through 1D and 2D infrared spectroscopy, especially using isotope labeling, as already demonstrated in the literature for peptides of similar length,89−93 not only for the strength of ion-binding but also for the peptoid parameters.





(1) Kunz, W. Specific Ion Effects in Colloidal and Biological Systems. Curr. Opin. Colloid Interface Sci. 2010, 15, 34−39. (2) Tobias, D. J.; Stern, A. C.; Baer, M. D.; Levin, Y.; Mundy, C. J. Simulation and Theory of Ions at Atmospherically Relevant Aqueous Liquid-Air Interfaces. Annu. Rev. Phys. Chem. 2013, 64, 339−359. (3) Katz, A. K.; Glusker, J. P.; Beebe, S. A.; Bock, C. W. Calcium Ion Coordination: A Comparison with that of Beryllium, Magnesium, and Zinc. J. Am. Chem. Soc. 1996, 118, 5752−5763. (4) Okur, H. I.; Kherb, J.; Cremer, P. S. Cations Bind Only Weakly to Amides in Aqueous Solutions. J. Am. Chem. Soc. 2013, 135, 5062− 5067. (5) Pluharova, E.; Baer, M. D.; Mundy, C. J.; Schmidt, B.; Jungwirth, P. Aqueous Cation-Amide Binding: Free Energies and IR Spectral Signatures by Ab Initio Molecular Dynamics. J. Phys. Chem. Lett. 2014, 5, 2235−2240. (6) Bertrand, M.; Sy, D.; Brack, A. Salt-dependent Conformational Diversity of Alternating Poly(Glu-Leu). J. Pept. Res. 1997, 49, 269− 272. (7) Hollosi, M.; Otvos, L.; Urge, L.; Kajtar, J.; Perczel, A.; Laczko, I.; Vadasz, Z. S.; Fasman, G. D. CA2+-Induced Conformational Transitions of Phosphorylated Peptides. Biopolymers 1993, 33, 497− 510. (8) Zarzhitsky, S.; Edri, H.; Azoulay, Z.; Cohen, I.; Ventura, Y.; Gitelman, A.; Rapaport, H. The Effect of pH and Calcium Ions on the Stability of Amphiphilic and Anionic beta-Sheet Peptide Hydrogels. Biopolymers 2013, 100, 760−772. (9) Sorensen, A. B.; Sondergaard, M. T.; Overgaard, M. T. Calmodulin in a Heartbeat. FEBS J. 2013, 280, 5511−5532. (10) Levin, Y.; dos Santos, A. P.; Diehl, A. Ions at the Air-Water Interface: An End to a Hundred-Year-Old Mystery? Phys. Rev. Lett. 2009, 103, 257802. (11) Levin, Y. Polarizable Ions at Interfaces. Phys. Rev. Lett. 2009, 102, 084715. (12) Lewis, T.; Winter, B.; Stern, A. C.; Baer, M. D.; Mundy, C. J.; Tobias, D. J.; Hemminger, J. C. Does Nitric Acid Dissociate at the Aqueous Solution Surface? J. Phys. Chem. C 2011, 115, 21183−21190. (13) Baer, M. D.; Stern, A. C.; Levin, Y.; Tobias, D. J.; Mundy, C. J. Electrochemical Surface Potential Due to Classical Point Charge Models Drives Anion Adsorption to the Air-Water Interface. J. Phys. Chem. Lett. 2012, 3, 1565−1570. (14) Stern, A. C.; Baer, M. D.; Mundy, C. J.; Tobias, D. J. Thermodynamics of Iodide Adsorption at the Instantaneous Air-water Interface. J. Chem. Phys. 2013, 138, 114709. (15) Baer, M. D.; Mundy, C. J. An Ab Initio Approach to Understanding the Specific Ion Effect. Faraday Discuss. 2013, 160, 89−101. (16) Baer, M. D.; Kuo, I.-F. W.; Tobias, D. J.; Mundy, C. J. Toward a Unified Picture of the Water Self-Ions at the Air-Water Interface: A Density Functional Theory Perspective. J. Phys. Chem. B 2014, 118, 8364−8372. (17) Levin, Y.; dos Santos, A. P. Ions at Hydrophobic Interfaces. J. Phys.: Condens. Matter 2014, 26, 203101. (18) Duignan, T. T.; Parsons, D. F.; Ninham, B. W. Ion Interactions with the Air-Water Interface Using a Continuum Solvent Model. J. Phys. Chem. B 2014, 118, 8700−8710. (19) dos Santos, A. P.; Levin, Y. Ions at the Water-oil Interface: Interfacial Tension of Electrolyte Solutions. Langmuir 2012, 28, 1304− 1308.

ASSOCIATED CONTENT

S Supporting Information *

This material is available free of charge via the Internet at http://pubs.acs.org/. The Supporting Information is available free of charge on the ACS Publications website at DOI: 10.1021/acs.jpcb.5b12277. A zip-file containing topology, structure, and input files for the classical FF simulations. (ZIP)



REFERENCES

AUTHOR INFORMATION

Corresponding Author

*E-mail: [email protected]. Notes

The authors declare no competing financial interest.



ACKNOWLEDGMENTS M.D.D. was supported by MS3 (Materials Synthesis and Simulation Across Scales) Initiative at Pacific Northwest National Laboratory. Research was funded by the Laboratory Directed Research and Development program at Pacific Northwest National Laboratory. PNNL is a multiprogram 2206

DOI: 10.1021/acs.jpcb.5b12277 J. Phys. Chem. B 2016, 120, 2198−2208

Article

The Journal of Physical Chemistry B (20) dos Santos, A. P.; Levin, Y. Surface and Interfacial Tensions of Hofmeister Electrolytes. Faraday Discuss. 2013, 160, 75−87. (21) Otten, D. E.; Shaffer, P. R.; Geissler, P. L.; Saykally, R. J. Elucidating the Mechanism of Selective Ion Adsorption to the Liquid Water Surface. Proc. Natl. Acad. Sci. U. S. A. 2012, 109, 701−705. (22) Saxena, A.; Garcia, A. E. Multisite Ion Model in Concentrated Solutions of D ivalent Cations (MgCl2 and CaCl2): Osmotic Pressure Calculations. J. Phys. Chem. B 2015, 119, 219−227. (23) Rousseau, R.; Schenter, G. K.; Fulton, J. L.; Linehan, J. C.; Engelhard, M. H.; Autrey, T. Defining Active Catalyst Structure and Reaction Pathways from ab Initio Molecular Dynamics and Operando XAFS: Dehydrogenation of Dimethylaminoborane by Rhodium Clusters. J. Am. Chem. Soc. 2009, 131, 10516−10524. (24) Fulton, J. L.; Schenter, G. K.; Baer, M. D.; Mundy, C. J.; Dang, L. X.; Balasubramanian, M. Probing the Hydration Structure of Polarizable Halides: A Multiedge XAFS and Molecular Dynamics Study of the Iodide Anion. J. Phys. Chem. B 2010, 114, 12926−12937. (25) Cauet, E.; Bogatko, S.; Weare, J. H.; Fulton, J. L.; Schenter, G. K.; Bylaska, E. J. Structure and Dynamics of the Hydration Shells of the Zn2+ Ion from Ab Initio Molecular Dynamics and Combined Ab Initio and Classical Molecular Dynamics Simulations. J. Chem. Phys. 2010, 132, 194502. (26) D’Angelo, P.; Migliorati, V.; Guidoni, L. Hydration Properties of the Bromide Aqua Ion: the Interplay of First Principle and Classical Molecular Dynamics, and X-ray Absorption Spectroscopy. Inorg. Chem. 2010, 49, 4224−4231. (27) Baer, M. D.; Pham, V.-T.; Fulton, J. L.; Schenter, G. K.; Balasubramanian, M.; Mundy, C. J. Is Iodate a Strongly Hydrated Cation? J. Phys. Chem. Lett. 2011, 2, 2650−2654. (28) Fulton, J. L.; Bylaska, E. J.; Bogatko, S.; Balasubramanian, M.; Cauet, E.; Schenter, G. K.; Weare, J. H. Near-Quantitative Agreement of Model-Free DFT-MD Predictions with XAFS Observations of the Hydration Structure of Highly Charged Transition-Metal Ions. J. Phys. Chem. Lett. 2012, 3, 2588−2593. (29) Bogatko, S.; Cauet, E.; Bylaska, E.; Schenter, G.; Fulton, J.; Weare, J. The Aqueous Ca2+ System, in Comparison with Zn2+, Fe3+, and Al3+: An Ab Initio Molecular Dynamics Study. Chem. - Eur. J. 2013, 19, 3047−3060. (30) Baer, M. D.; Fulton, J. L.; Balasubramanian, M.; Schenter, G. K.; Mundy, C. J. Persistent Ion Pairing in Aqueous Hydrochloric Acid. J. Phys. Chem. B 2014, 118, 7211−7220. (31) Sessa, F.; D’Angelo, P.; Guidoni, L.; Migliorati, V. Hidden Hydration Structure of Halide Ions: an Insight into the Importance of Lone Pairs. J. Phys. Chem. B 2015, 119, 15729−15737. (32) Timko, J.; De Castro, A.; Kuyucak, S. Ab Initio Calculation of the Potential of Mean Force for Dissociation of Aqueous Ca-Cl. J. Chem. Phys. 2011, 134, 204510. (33) Pluharova, E.; Marsalek, O.; Schmidt, B.; Jungwirth, P. Ab Initio Molecular Dynamics Approach to a Quantitative Description of Ion Pairing in Water. J. Phys. Chem. Lett. 2013, 4, 4177−4181. (34) Allolio, C.; Salas-Illanes, N.; Desmukh, Y. S.; Hansen, M. R.; Sebastiani, D. H-Bonding Competition and Clustering in Aqueous LiI. J. Phys. Chem. B 2013, 117, 9939−9946. (35) Kahlen, J.; Salimi, L.; Sulpizi, M.; Peter, C.; Donadio, D. Interaction of Charged Amino-Acid Side Chains with Ions: An Optimization Strategy for Classical Force Fields. J. Phys. Chem. B 2014, 118, 3960−3972. (36) Pluharova, E.; Baer, M. D.; Schenter, G. K.; Jungwirth, P.; Mundy, C. J. The Dependence of the Rate of LiF Ion Pairing on the Description of Molecular Interaction. J. Phys. Chem. B 2015, DOI: 10.1021/acs.jpcb.5b09344. (37) Pavlov, M.; Siegbahn, P.; Sandstrom, M. Hydration of Beryllium, Magnesium, Calcium, and Zinc Ions Using Density Functional Theory. J. Phys. Chem. A 1998, 102, 219−228. (38) Leontyev, I.; Stuchebrukhov, A. Accounting for Electronic Polarization in Non-Polarizable Force Fields. Phys. Chem. Chem. Phys. 2011, 13, 2613−2626.

(39) Kohagen, M.; Mason, P. E.; Jungwirth, P. Accurate Description of Calcium Solvation in Concentrated Aqueous Solutions. J. Phys. Chem. B 2014, 118, 7902−7909. (40) Kohagen, M.; Pluharova, E.; Mason, P. E.; Jungwirth, P. Exploring Ion-Ion Interactions in Aqueous Solutions by a Combination of Molecular Dynamics and Neutron Scattering. J. Phys. Chem. Lett. 2015, 6, 1563−1567. (41) Kohagen, M.; Lepsik, M.; Jungwirth, P. Calcium Binding to Calmodulin by Molecular Dynamics with Effective Polarization. J. Phys. Chem. Lett. 2014, 5, 3964−3969. (42) Yoo, B.; Kirshenbaum, K. Peptoid Architectures: Elaboration, Actuation, and Application. Curr. Opin. Chem. Biol. 2008, 12, 714− 721. (43) Laursen, J. S.; Harris, P.; Fristrup, P.; Olsen, C. A. Triangular Prism-shaped Beta-peptoid Helices as Unique Biomimetic Scaffolds. Nat. Commun. 2015, 6, 7013. (44) Laursen, J. S.; Engel-Andreasen, J.; Olsen, C. A. beta-Peptoid Foldamers at Last. Acc. Chem. Res. 2015, 48, 2696−2704. (45) Nam, K. T.; Shelby, S. A.; Choi, P. H.; Marciel, A. B.; Chen, R.; Tan, L.; Chu, T. K.; Mesch, R. A.; Lee, B. C.; Connolly, M. D.; et al. Free-floating Ultrathin Two-dimensional Crystals from Sequencespecific Peptoid Polymers. Nat. Mater. 2010, 9, 454−60. (46) Mannige, R. V.; Haxton, T. K.; Proulx, C.; Robertson, E. J.; Battigelli, A.; Butterfoss, G. L.; Zuckermann, R. N.; Whitelam, S. Peptoid Nanosheets Exhibit a new Secondary-Structure Motif. Nature 2015, 526, 415−+. (47) Kirshenbaum, K.; Barron, A. E.; Goldsmith, R. A.; Armand, P.; Bradley, E. K.; Truong, K. T. V.; Dill, K. A.; Cohen, F. E.; Zuckermann, R. N. Sequence-specific Polypeptoids: A Diverse Family of Heteropolymers with Stable Secondary Structure. Proc. Natl. Acad. Sci. U. S. A. 1998, 95, 4303−4308. (48) Sanborn, T. J.; Wu, C. W.; Zuckermann, R. N.; Barron, A. E. Extreme Stability of Helices Formed by Water-soluble Poly-Nsubstituted Glycines (polypeptoids) with alpha-chiral Side Chains. Biopolymers 2002, 63, 12−20. (49) Sun, J.; Zuckermann, R. N. Peptoid Polymers: A Highly Designable Bioinspired Material. ACS Nano 2013, 7, 4715−4732. (50) Culf, A. S.; Ouellette, R. J. Solid-Phase Synthesis of NSubstituted Glycine Oligomers (alpha-Peptoids) and Derivatives. Molecules 2010, 15, 5282−5335. (51) Zuckermann, R. N.; Kodadek, T. Peptoids as Potential Therapeutics. Current Opinion in Molecular Therapeutics 2009, 11, 299−307. (52) Goodman, C. M.; Choi, S.; Shandler, S.; DeGrado, W. F. Foldamers as Versatile Frameworks for the Design and Evolution of Function. Nat. Chem. Biol. 2007, 3, 252−262. (53) Matsuurua, K. Rational Design of Self-assembled Proteins and Peptides for Nano- and Micro-sized Architectures. RSC Adv. 2014, 4, 2942−2953. (54) Tymoczko, J. L.; Berg, J. M.; Stryer, L. Biochemistry. Second printing; WH Freeman and Company: New York, 2010. (55) Wang, J. M.; Wolf, R. M.; Caldwell, J. W.; Kollman, P. A.; Case, D. A. Development and Testing of a General AMBER Force Field. J. Comput. Chem. 2004, 25, 1157−1174. (56) Wang, J.; Wang, W.; Kollman, P. A.; Case, D. A. Automatic Atom Type and Bond Type Perception in Molecular Mechanical Calculations. J. Mol. Graphics Modell. 2006, 25, 247−260. (57) Jakalian, A.; Bush, B. L.; Jack, D. B.; Bayly, C. I. Fast, Efficient Generation of High-quality Atomic Charges. AM1-BCC Model: I. Method J. Comput. Chem. 200021132−146.10.1002/(SICI)1096987X(20000130)21:23.3.CO;2-G (58) Jakalian, A.; Jack, D. B.; Bayly, C. I. Fast, Efficient Generation of High-quality Atomic Charges. AM1-BCC Model: II. Parameterization and Validation. J. Comput. Chem. 2002, 23, 1623−1641. (59) Voelz, V. A.; Dill, K. A.; Chorny, I. Peptoid Conformational Free Energy landscapes From Implicit-Solvent Molecular Simulations in AMBER. Biopolymers 2011, 96, 639−650. (60) Butterfoss, G. L.; Yoo, B.; Jaworski, J. N.; Chorny, I.; Dill, K. A.; Zuckermann, R. N.; Bonneau, R.; Kirshenbaum, K.; Voelz, V. A. De 2207

DOI: 10.1021/acs.jpcb.5b12277 J. Phys. Chem. B 2016, 120, 2198−2208

Article

The Journal of Physical Chemistry B Novo Structure Prediction and Experimental Characterization of Folded Peptoid Oligomers. Proc. Natl. Acad. Sci. U. S. A. 2012, 109, 14320−14325. (61) Mukherjee, S.; Zhou, G.; Michel, C.; Voelz, V. A. Insights into Peptoid Helix Folding Cooperativity from an Improved Backbone Potential. J. Phys. Chem. B 2015, 119, 15407−15417. (62) Mirijanian, D. T.; Mannige, R. V.; Zuckermann, R. N.; Whitelam, S. Development and Use of an Atomistic CHARMMBased Forcefield for Peptoid Simulation. J. Comput. Chem. 2014, 35, 360−370. (63) Sousa da Silva, A. W.; Vranken, W. F. ACPYPE - AnteChamber PYthon Parser interfacE. BMC Res. Notes 2012, 5, 367−367. (64) Duan, Y.; Wu, C.; Chowdhury, S.; Lee, M. C.; Xiong, G. M.; Zhang, W.; Yang, R.; Cieplak, P.; Luo, R.; Lee, T.; et al. A Point-charge Force Field for Molecular Mechanics Simulations of Proteins based on Condensed-phase Quantum Mechanical Calculations. J. Comput. Chem. 2003, 24, 1999−2012. (65) Auffinger, P.; Cheatham, T. E.; Vaiana, A. C. Spontaneous Formation of KCl Aggregates in Biomolecular Simulations: A Force Field Issue? J. Chem. Theory Comput. 2007, 3, 1851−1859. (66) Dang, L. X.; Kollman, P. A. Free Energy of Association of the K +:18-Crown-6 Complex in Water: A New Molecular Dynamics Study. J. Phys. Chem. 1995, 99, 55−58. (67) Berendsen, H. J. C.; Postma, J. P. M.; van Gunsteren, W. F.; Hermans, J. In In Intermolecular Forces; Pullman, B., Ed.; Reidel: Dordrecht, 1981. (68) Onufriev, A.; Bashford, D.; Case, D. A. Exploring Protein Native States and Large-scale Conformational Changes with a Modified Generalized Born Model. Proteins: Struct., Funct., Genet. 2004, 55, 383−394. (69) Bussi, G.; Donadio, D.; Parrinello, M. Canonical Sampling Through Velocity Rescaling. J. Chem. Phys. 2007, 126, 014101. (70) Pronk, S.; Pall, S.; Schulz, R.; Larsson, P.; Bjelkmar, P.; Apostolov, R.; Shirts, M. R.; Smith, J. C.; Kasson, P. M.; van der Spoel, D.; et al. GROMACS 4.5: a High-throughput and Highly Parallel Open Source Molecular Simulation Toolkit. Bioinformatics 2013, 29, 845−854. (71) Darden, T.; York, D.; Pedersen, L. Particle Mesh Ewald - An N.log(n) Method For Ewald Sums In Large Systems. J. Chem. Phys. 1993, 98, 10089−10092. (72) Parrinello, M.; Rahman, A. Polymorphic Transitions In Singlecrystals - A New Molecular-Dynamics Method. J. Appl. Phys. 1981, 52, 7182−7190. (73) Hoover, W. G. Canonical Dynamics - Equilibrium Phase-space Distributions. Phys. Rev. A: At., Mol., Opt. Phys. 1985, 31, 1695−1697. (74) Grossfield, A. WHAM: the Weighted Histogram Analysis Method. 2011; http://membrane.urmc.rochester.edu/content/wham. (75) Hartigan, J. A.; Wong, M. A. Algorithm AS 136: A K-Means Clustering Algorithm. Journal of the Royal Statistical Society. Series C (Applied Statistics) 1979, 28, 100−108. (76) R Core Team, R: A Language and Environment for Statistical Computing; R Foundation for Statistical Computing: Vienna, Austria, 2013. (77) VandeVondele, J.; Krack, M.; Mohamed, F.; Parrinello, M.; Chassaing, T.; Hutter, J. QUICKSTEP: Fast and Accurate Density Functional Calculations using a Mixed Gaussian and Plane Waves Approach. Comput. Phys. Commun. 2005, 167, 103−128. (78) VandeVondele, J.; Hutter, J. Gaussian Basis Sets for Accurate Calculations on Molecular Systems in Gas and Condensed Phases. J. Chem. Phys. 2007, 127, 114105. (79) Goedecker, S.; Teter, M.; Hutter, J. Separable Dual-Space Gaussian Pseudopotentials. Phys. Rev. B: Condens. Matter Mater. Phys. 1996, 54, 1703−1710. (80) Martyna, G. J.; Klein, M. L.; Tuckerman, M. E. Nose-Hoover Chains - The Canonical Ensemble via Continous Dynamics. J. Chem. Phys. 1992, 97, 2635−2643. (81) Becke, A. D. Density-Functional Exchange-Energy Approximation with Correct Asymptotic Behavior. Phys. Rev. A: At., Mol., Opt. Phys. 1988, 38, 3098−3100.

(82) Lee, C.; Yang, W.; Parr, R. G. Development of the Colle-Salvetti Correlation-Energy Formula into a Functional of the Electron Density. Phys. Rev. B: Condens. Matter Mater. Phys. 1988, 37, 785−789. (83) Grimme, S. Accurate Description of Van der Waals Complexes by Density Functional Theory Including Empirical Corrections. J. Comput. Chem. 2004, 25, 1463−1473. (84) Dang, L. X.; Truong, T. B.; Ginovska-Pangovska, B. Note: Interionic Potentials of Mean Force for Ca2+-Cl- in Polarizable Water. J. Chem. Phys. 2012, 136, 126101. (85) Dudev, T.; Lim, C. J. Phys. Chem. B 2004, 108, 4546−4557. (86) Pham, V.-T.; Fulton, J. L. Ion-pairing in Aqueous CaCl2 and RbBr Solutions: Simultaneous Structural Refinement of XAFS and XRD Data. J. Chem. Phys. 2013, 138, 044201. (87) Baer, M. D.; Mundy, C. J. Local Aqueous Solvation Structure Around Ca2+ During Ca2+−Cl− Pair Formation. J. Phys. Chem. B 2016, DOI: 10.1021/acs.jpcb.5b09579. (88) Ikeda, T.; Boero, M.; Terakura, K. Hydration Properties of Magnesium and Calcium Ions from Constrained First Principles Molecular Dynamics. J. Chem. Phys. 2007, 127, 074503. (89) Brewer, S. H.; Song, B.; Raleigh, D. P.; Dyer, R. B. Residue Specific Resolution of Protein Folding Dynamics using Isotope-edited Infrared Temperature Jump Spectroscopy. Biochemistry 2007, 46, 3279−3285. (90) Hauser, K.; Krejtschi, C.; Huang, R.; Wu, L.; Keiderling, T. A. Site-specific Relaxation Kinetics of a Tryptophan Zipper Hairpin Peptide using Temperature-jump IR Spectroscopy and Isotopic Labeling. J. Am. Chem. Soc. 2008, 130, 2984−2992. (91) Shim, S.-H.; Gupta, R.; Ling, Y. L.; Strasfeld, D. B.; Raleigh, D. P.; Zanni, M. T. Two-dimensional IR spectroscopy and Isotope Labeling Defines the Pathway of Amyloid Formation with Residuespecific Resolution. Proc. Natl. Acad. Sci. U. S. A. 2009, 106, 6614− 6619. (92) Decatur, S. M. IR Spectroscopy of Isotope-labeled Helical Peptides: Probing the Effect of N-acetylation on Helix Stability. Biopolymers 2000, 54, 180−185. (93) Lessing, J.; Roy, S.; Reppert, M.; Baer, M.; Marx, D.; Jansen, T. L. C.; Knoester, J.; Tokmakoff, A. Identifying Residual Structure in Intrinsically Disordered Systems: A 2D IR Spectroscopic Study of the GVGXPGVG Peptide. J. Am. Chem. Soc. 2012, 134, 5032−5035.

2208

DOI: 10.1021/acs.jpcb.5b12277 J. Phys. Chem. B 2016, 120, 2198−2208