Solvation of the Guanidinium Ion in Pure Aqueous Environments: A

Subscriber access provided by University of Florida | Smathers Libraries

Article

Solvation of the Guanidinium Ion in Pure Aqueous Environments. A Theoretical Study from an "Ab Initio"-based Polarizable Force Field Céline Houriez, Michael Mautner, and Michel Masella J. Phys. Chem. B, Just Accepted Manuscript • DOI: 10.1021/acs.jpcb.7b07874 • Publication Date (Web): 28 Nov 2017 Downloaded from http://pubs.acs.org on December 6, 2017

Just Accepted “Just Accepted” manuscripts have been peer-reviewed and accepted for publication. They are posted online prior to technical editing, formatting for publication and author proofing. The American Chemical Society provides “Just Accepted” as a free service to the research community to expedite the dissemination of scientific material as soon as possible after acceptance. “Just Accepted” manuscripts appear in full in PDF format accompanied by an HTML abstract. “Just Accepted” manuscripts have been fully peer reviewed, but should not be considered the official version of record. They are accessible to all readers and citable by the Digital Object Identifier (DOI®). “Just Accepted” is an optional service offered to authors. Therefore, the “Just Accepted” Web site may not include all articles that will be published in the journal. After a manuscript is technically edited and formatted, it will be removed from the “Just Accepted” Web site and published as an ASAP article. Note that technical editing may introduce minor changes to the manuscript text and/or graphics which could affect content, and all legal disclaimers and ethical guidelines that apply to the journal pertain. ACS cannot be held responsible for errors or consequences arising from the use of information contained in these “Just Accepted” manuscripts.

The Journal of Physical Chemistry B is published by the American Chemical Society. 1155 Sixteenth Street N.W., Washington, DC 20036 Published by American Chemical Society. Copyright © American Chemical Society. However, no copyright claim is made to original U.S. Government works, or works produced by employees of any Commonwealth realm Crown government in the course of their duties.

Page 1 of 37 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

Solvation of the Guanidinium Ion in Pure Aqueous Environments. A Theoretical Study from an “Ab Initio”-Based Polarizable Force Field. Céline Houriez,† Michael Meot-Ner (Mautner),‡ and Michel Masella∗,¶ †MINES ParisTech, PSL Research University, CTP - Centre Thermodynamique des Procédés, 35 rue Saint-Honoré, 77300 Fontainebleau, France ‡Department of Chemistry, Virginia Commonwealth University, Richmond, Virginia 23284-2006, United States, and Department of Chemistry, University of Canterbury, Christchurch, New Zealand 8001 ¶Laboratoire de Biologie Structurale et Radiobiologie, Service de Bioénergétique, Biologie Structurale et Mécanismes, Institut de Biologie et de Technologies de Saclay, CEA Saclay, F-91191 Gif sur Yvette Cedex, France E-mail: [email protected]

1

ACS Paragon Plus Environment

The Journal of Physical Chemistry 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Abstract

1

2

We report simulation results regarding the hydration process of the guanidinium

3

cation in water droplets and in bulk liquid water, at a low concentration of 0.03 M,

4

performed using a polarizable approach to model both water/water and ion/water

5

interactions. In line with earlier theoretical studies, our simulations show a preferential

6

orientation of guanidinium at water/vacuum interfaces, i.e. a parallel orientation of the

7

guanidinium plane to the aqueous surface. In an apparent contradiction with earlier

8

simulation studies, we show also that guanidinium has a stronger propensity for the

9

cores of aqueous systems than the ammonium cation. However our bulk simulation

10

conditions correspond to weaker cation concentrations than in earlier studies, by two

11

orders of magnitude, and that the same simulations performed using a standard non-

12

polarizable force field leads to the same conclusion. From droplet data, we extrapolate

13

the guanidinium single hydration enthalpy value to be -82.9 ± 2.2 kcal mol−1 . That

14

is about half as large as the sole experimental estimate reported to date, about -

15

144 kcal mol−1 . Our result yields a guanidinium absolute bulk hydration free energy

16

at ambiant conditions to be -78.4 ± 2.6 kcal mol−1 , a value smaller by 3 kcal mol−1

17

compared to ammonium. The relatively large magnitude of our guanidinium hydration

18

free energy estimate suggests the Gdm+ protein denaturing properties to result from a

19

competition between the cation hydration effects and the cation/protein interactions, a

20

competition that can be modulated by weak differences in the protein or in the cation

21

chemical environment.

2


Page 2 of 37

Page 3 of 37 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


22

1

Introduction

23

+ The guanidinium ion C(NH2 )+ 3 , denoted by Gdm , is an overall large, symmetric and planar

24

organic cation that plays a key role in biochemistry, in particular because it is the main

25

constituent of the arginine amino acid side chain. Gdm+ salts are also well known protein

26

denaturing agents, 1 even if we may note the conflicting observation that Gdm+ can be found

27

both at the head and at the tail of the cation Hofmeister series, 2,3 i.e. a series where ions are

28

ordered based on their ability in stabilizing/destabilizing folded protein structures. Other

29

striking properties of Gdm+ are (1) its ability to form charge-like pairs in aqueous environ-

30

ments, 3–5 even when belonging to arginine side chains, 6,7 and (2) its ability to preferentially

31

interact with aromatics than with aliphatics in aqueous phase. 8 The original properties of

32

Gdm+ are inferred to arise from asymmetric hydration properties tied to its planar structure.

33

For instance, the three Gmd+ NH2 moieties defining the cation symmetry plane are shown

34

from neutron scattering experiments to form NH − O hydrogen bonds with water molecules

35

in the Gdm+ first hydration shell, whereas the region above the Gdm+ symmetry plane is

36

considered as hydrophobic. 4

37

Despite noticeable efforts devoted to understand the original properties of Gdm+ , the

38

role of Gdm+ /water interactions on these properties is still largely debated, in particular

39

because of the lack of reliable experimental evidence and data regarding the Gdm+ hydration

40

process. For instance, the sole available experimental-based estimate of the Gdm+ single ion

41

hydration enthalpy ∆Hhyd (Gdm+ ) is 144 kcal mol−1 as reported recently by Y. Marcus. 9

42

However that value is considered as “puzzling” by the latter author 9 as it is even larger

43

by 10 kcal mol−1 than all the reported hydration enthalpy estimates of the smallest alkali

44

cation Li+ and also because that largely conflicts with the well accepted property of Gdm+ ,

45

i.e. it interacts weakly with water. 10–12 We may also note here that the latter experimental

46

∆Hhyd (Gdm+ ) value is conflicting with theoretical hydration free energy estimates ∆Ghyd

47

for Gdm+ that are about 60 kcal mol−1 from quantum computations where the solvent is

48

accounted for using a polarizable continuum approach. 13 Assuming the reliability of the 3



49

latter experimental and theoretical thermodynamic values, they yield a hydration entropy

50

component −T∆Shyd for Gdm+ of 84 kcal mol−1 , a value which is one order of magnitude

51

larger than all the reported values for alkali and small organic cations 14 as well as the

52

experimental-based value reported by Marcus 9 for Gdm+ (about 5 ± 1 kcal mol−1 ). Lastly

53

we may here also quote a recent experimental study of Heiles et al 15 dealing with the solvation

54

of Gdm+ in small water droplets (comprising at most 100 water molecules) that shows the

55

strength of the Gdm+ /water interactions to be modulated by the chemical surroundings.

56

In particular the latter authors concluded that Gdm+ interacts more strongly with water

57

in chemical environments where water is partially excluded and that this feature is tied

58

to the effectiveness of Gdm+ as protein denaturant. Nevertheless these findings is a priori

59

contradictory with the large magnitude of the ∆Hhyd (Gdm+ ) value in bulk water.

60

Theoretically, most of the studies devoted to investigating the Gdm+ hydration process

61

were performed using standard pairwise force fields, 4,5,7,16–19 whose accuracy is still debated

62

for such a cation. 17,18 More sophisticated theoretical approaches have been proposed, based

63

on a quantum ab initio scheme 20 or on polarizable force fields, 21–23 for instance. However,

64

these approaches have still not been used to perform an in-depth analysis of the Gdm+

65

hydration process in various aqueous environments. Recently we proposed a rigid polarizable

66

force field for water, TCPE/2013, 24 whose parameters are mostly assigned from high level

67

quantum computations regarding small water clusters. In particular, this model is a priori

68

well suited to be used for investigating ion hydration processes in aqueous environments

69

(small clusters, large droplets, water/vacuum interfaces and bulk phase) not only because it

70

is able to accurately model a wide range of water systems under different physical conditions,

71

but also because we paid attention to assign its parameters to reproduce the magnitude of

72

the repulsive water/water interactions in cation first hydration shells. We have already

73

applied this water model to investigate the hydration process in water finite size systems or

74

in bulk water of monoatomic heavy cations (Cm3+ and Th4+ ), 24 of methylated ammonium

75

cations, 25,26 of monoatomic anions 27 and of alkylated carboxylates. 28 In these studies, we

4


Page 4 of 37

Page 5 of 37 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


76

only assigned the force field parameters handling the ion/water interactions on the basis

77

of high level quantum computational results. Even if there is still some room to improve

78

our force-field approach, most of the results are in line with experiment. For instance,

79

by extrapolating ion/water interaction results from nanodrop simulations, we computed an

80

estimate of the proton hydration enthalpy ∆Hhyd (H+ ) in liquid water that ranges from 269

81

to 274 kcal mol−1 , depending on the ion considered, in good agreement with experimental

82

based estimates (ranging from 271 to 275 kcal mol−1 29–31 ). This demonstrates the capacity

83

of the so-called ab initio-based force-fields to provide accurate data regarding the hydration

84

process of any kind of ion, in particular those suffering from a lack of experimental evidence

85

in particular chemical environments.

86

The aim of the present paper is to apply our ab initio-based force field methodology to

87

the particular case of the Gdm+ ion to complement the available data regarding its hydration

88

process both in aqueous nanodrops and in bulk environments. At the difference of recent

89

theoretical studies dealing with guanidinium/chloride salts, 17,18 here we simulate Gdm+ as

90

solvated alone in water environnements in order to estimate its individual hydration prop-

91

erties. This is the first step to understand the hydration properties of aqueous Gdm+ /anion

92

salt solutions (in particular to discuss the existence and thus the possible consequences of

93

Gdm+ /anion/water cooperative effects).

94

As standard (i.e. non polarizable) force fields were shown to reproduce recent experi-

95

mental studies about the solvation of guanidinium/chloride salt solutions, 17,18 we performed

96

all our simulations twice by using also the standard force field proposed by Wernersson et

97

all 17 to cross check the reliability both of our polarizable approach and of standard force

98

fields to investigate the Gdm+ hydration properties. Lastly, as the ammonium cation is the

99

second kind of positively charged chemical groups present in amino acids, we will here also

100

compare the hydration properties of Gdm+ to the NH+ 4 ones as computed from out polariz-

101

able approach (here we consider the data reported in our recent studies 25,26 ) and also from

102

the standard additive force field used by Werner and coworkers 19 (however only in the liquid

5



103

water-vacuum case).

104

2

105

In the following, M is the total number of atoms considered, Mion is the total number of

106

atoms in an ion, Mµ is the total number of polarizable atoms, N is the total number of water

107

molecules and Mw is the total number of atoms belonging to the water molecules.

108

2.1

109

The interactions of Gdm+ with 1 to 3 water molecules are investigated using ab initio quan-

110

tum computations at the MP2 level (using the frozen core approximation) with the basis

111

sets aug-cc-pVXZ (X = D,T and Q). The Gdm+ /water binding energies, ∆U , are computed

112

from aggregate geometries optimized at the same level of theory for X=D and T, and from

113

aug-cc-pVTZ ones for X=Q. All ab initio computations are performed with the package of

114

programs GAUSSIAN09 version C01. 32 Complete Basis Set (CBS) limit estimates of the

115

Gdm+ /water aggregate binding energies are extrapolated from our MP2/aug-c–pVXZ data

116

according to standard schemes. 33,34 We also computed at the MP2/CBS limit level of the-

117

ory two radial elongation profiles corresponding to the interaction of Gdm+ with a single

118

water molecule. These profiles correspond to (1) a water molecule interacting with two NH2

119

moieties in the Gdm+ plane (see the cation/water hetero dimer structure shown in Figure

120

1), and (2) a water molecule set orthogonally to the Gdm+ plane (see Figure 2). These

121

two profiles are computed for oxygen-carbon distances ranging from 2 to 8 ˚ A and regularly

122

spaced. We also performed single energy point calculations at the coupled cluster level of

123

theory, including single, double and triple excitations, CCSD(T), using the above basis sets

124

for a few hydrated Gdm+ aggregates. The CCSD(T) results agree with the MP2 ones within

125

less than 0.1 kcal mol−1 .

Theoretical details

Quantum Computations

6


Page 6 of 37

Page 7 of 37 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


126

2.2

Molecular Force Fields

127

Our polarizable force field is based on a decomposition of the total potential energy U of

128

Gdm+ /water systems into three main energy terms, namely the ion internal energy U rel ,

129

and the ion/water U iw and water/water U ww interaction energies. For U ww , we consider

130

the rigid and polarizable water model TCPE/2013 24 that was shown to provide an accurate

131

description of both isolated water clusters in vacuum and bulk water systems for a wide range

132

of temperature and pressure conditions. Besides a specific anisotropic many body energy

133

term U lh introduced to accurately model water hydrogen bond (HB) networks, TCPE/2013

134

also considers a repulsive U rep term, a Coulombic U qq term, and a polarization U pol term

135

based on an induced dipole moment approach that includes damping effects: 24

0

U rep =

X

aij exp (−bij rij ) ,

(1)

1 X qi qj , 4π0 i,j rij

(2)

i,j 136 0

U qq = 137

N

U pol

N

N

N∗

µ µ µ µ X p2i X 1X 1X pi Tij pj . = − pi · Eqi − 2 i=1 αi i=1 2 i=1 j=1

(3)

138

Here, the sums run on pairs of atoms i and j belonging to different molecules. rij is the

139

distance between the atoms i and j. The {qi } are the static charges located on atomic centers.

140

aij and bij are adjustable parameters. αi is the isotropic polarizability of a polarizable atomic

141

centre i (all the atoms except the hydrogens), Tij is the dipolar interaction tensor and Eqi

142

is the electric field generated on the polarizable center i by the surrounding static charges

143

qj . Tij and Eqi include both a short range damping component. These energy terms were

144

already detailed in our former studies dealing with ammonium and carboxylate ions. 25,28 To

145

model the Gdm+ /water interactions, we consider also the three above intermolecular terms

7



146

Page 8 of 37

to which we add a standard additive dispersion term U disp

U

disp

=−

X i,j

ij

∗ rij rij

!6

.

(4)

147

∗ are two adjustable parameters. The Gdm+ /water interaction potential U iw Here, ij and rij

148

is thus the same sum of four terms that we considered to investigate the hydration process

149

of methylated ammonium cations 0

U iw = U rep + U qq + U pol + U disp .

(5)

150

The Gdm+ intramolecular degrees of freedom are handled by a dedicated relaxation en-

151

ergy term U rel , which includes standard stretching, bending and torsional terms and an

152

improper torsional term U imp =

153

gles centered on the cation carbon and nitrogen atoms. As in our former studies, we still

154

consider for the stretching and bending intramolecular term the parameters of the force

155

field CHARMM v2.7. 35 As we systematically constrain the X − H bonds and the H − X − H

156

angles along our simulations and as we simulate Gdm+ in pure water environments (that

157

prevents the force-field artifacts mentioned in Ref. 17 and originating from Gdm+ /chloride

158

interactions), the choice of the stretching and bending parameters has an overall negligible

159

effect on our gas phase cluster results. Concerning the standard dihedral 6 N − C − N − H

160

torsional term, we assign its parameters to reproduce the ab initio energy profile correspond-

161

ing to that degree of freedom from MP2/aug-cc-pVTZ data (see Supporting Information).

1 k ψ2, 2 imp

where ψ stands for the improper dihedral an-

162

To assign the U iw parameters, we consider as reference data the ab initio results con-

163

cerning three small guanidium/(H2 O)n clusters (n=1-3) in their geometries optimized at the

164

MP2/aug-cc-pVTZ level and the two radial elongation profiles mentioned in the above sec-

165

tion. All the energy data correspond to extrapolated MP2/CBS limit values from MP2(FC)/aug-

166

cc-pVXZ computations (X=D,T and Q). As in our former studies, 25,28 the static electrostatic

167

charges of the Gdm+ atoms were assigned from the quantum Natural Population Analysis 8


Page 9 of 37 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


168

results concerning the isolated ion at the MP2/aug-cc-pVTZ level. The Gdm+ atomic po-

169

larizabilities were assigned to reproduce the experimental ion molecular polarizability (4.4

170

˚ ˚3 ), A3 , 9 a value that agrees with the one computed at the MP2/aug-cc-pVTZ level, 4.0 A

171

˚3 , respectively. Here we by setting the carbon and nitrogen polarizabilities to 0.8 and 1.2 A

172

only consider the three Gdm+ nitrogens as dispersion centers. The parameter set is able to

173

reproduce the Gdm+ /water interaction energies within less than 0.5 kcal mol−1 on average

174

and the Gdm+ distances N · · · O within about 0.01 ˚ A (see Figure 1 and more details are

175

provided as Supporting Information).

176

For comparison purposes, we performed all our computations twice using a standard non

177

polarizable force field (for which the inter molecular interactions are decomposed into two

178

additive energy terms, the Lennard-Jones and the Coulombic ones) used in conjunction with

179

a parameter set that was shown to successfully reproduce experimental results regarding the

180

behavior of Gdm+ /Cl− salts at the liquid water-vacuum interface. 17,19 Regarding the stepwise

181

binding energies of the micro hydrated Gdm+ clusters shown Figure 1, the standard force

182

field provides fairly accurate results compared to quantum MP2/CBS data, within about 1

183

kcal mol−1 on average (see Figure 1). However it overestimates the Gdm+ nitrogen/water

184

oxygen distance in the latter clusters compared to quantum computations, by about 0.15 ˚ A

185

on average. Regarding the case of a water molecule interacting orthogonally to the Gdm+

186

plane, the standard force field predicts the carbon/oxygen distance corresponding to the

187

A. Gdm+ /water interaction energy minimum to be greater than the quantum one by 0.5 ˚

188

Note also that force field systematically overestimates the Gdm+ /water binding energy by

189

1 kcal mol−1 for carbon/oxygen distances greater than 3.5 ˚ A. Lastly, we investigated also

190

the behavior of the ammonium cation NH+ 4 at the liquid water/vacuum interface using the

191

standard additive force field used by the above authors. 17,36

9



192

2.3

Simulation details

193

We detail here briefly the theoretical protocol that we used to perform the simulations when

194

using both our polarizable force field and the standard one. More details are provided as

195

Supporting Information. All the molecular modeling computations were performed with

196

our in-house code POLARIS(MD). All the molecular dynamics, MD, simulations of ions

197

solvated either in nanodrops or at the liquid water-vacuum interface are performed in the

198

NVT ensemble, while bulk simulations are done in the NPT ensemble. The simulations

199

are all performed at the 10 ns scale according to the same theoretical protocol as in our

200

former studies. 25,26,28 In particular, bulk and liquid water-vacuum interface simulations are

201

performed using the SPME approach detailed in Ref. 37 The trajectories are sampled each 1

202

ps to generate the ensemble used to compute the statistical averages.

203

To minimize the effects of water evaporation phenomena in nanodrop simulations, we

204

consider the droplets confined in a spherical cavity, whose radius corresponds to the largest

205

distance droplet center of mass (COM)/water oxygen to which 12 ˚ A is added. As a water

206

molecule crosses the cavity boundary, it undergoes a smooth repulsive potential favoring

207

it returns to the main droplet region as we proposed recently. 26 Note that this repulsive

208

potential has no effect on the droplet dynamics, as it is zero for any atom/droplet COM

209

distance smaller than the cavity radius.

210

The potentials of mean force (PMF) corresponding to the interaction of Gdm+ with

211

a water system (droplet or liquid water) are computed according to a standard umbrella

212

sampling protocol. The degree of freedom harmonically constrained during the simulations

213

is the distance r between the Gdm+ carbon and the droplet COM or the projection z of the

214

distance between that carbon and the simulation cell center (SCC) on the axis orthogonal

215

to the liquid water-vacuum interface. We denote rc and zc the target values of the latter

216

parameters along the umbrella sampling simulations. The PMF are computed from the

217

constrained simulations at the 10 ns scale using the Umbrella Integration scheme. 38

10


Page 10 of 37

Page 11 of 37 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


218

3

Results and discussion

219

In the following, if not mentioned the results correspond to our polarizable force field. We

220

note l⊥ the axis orthogonal to the Gdm+ symmetry plane (computed as the cross product

221

of two of the three carbon-nitrogen vectors).

222

3.1

223

The Gdm+ carbon-water oxygen and Gdm+ nitrogen-water oxygen radial distribution func-

224

tions gC (r) and gN (r) in liquid water at ambient conditions (and the results of their inte-

225

gration) computed along our bulk MD simulations are shown in Figure 3. We also plot on

226

that figure the functions gX (r, φ0 ) corresponding to atom pairs for which the angle φ between

227

the axe connecting them and the vector l⊥ obeys |φ| ≥ φ0 , as well as the complementary

228

functions gX∗ (r, φ0 ) = gX (r) − gX (r, φ0 ). We set the angle φ0 to the values 30 and 60◦ for the

229

carbon and the nitrogen-centered functions, respectively. All these functions are scaled by

230

the mean water bulk density ρwater at ambiant conditions (0.0331 molecules per ˚ A3 ), and the

231

nitrogen ones are averaged over the three nitrogen data. The functions gC (r) and gN (r) con-

232

verge thus towards ρ˜water = 1 for large r values, and the functions gX∗ (r, φ0 ) converge towards

233

ρ˜φ0 = ρ˜water cos(φ0 ) (for φ0 = 30◦ and 60◦ we get then ρ˜φ0 = 0.87 and 0.50, respectively).

Guanidinium structural properties in bulk water

234

The function gC (r) is characterized by a first peak located at about 3.8 ˚ A that ends at

235

about 5.4 ˚ A. That peak corresponds to the first hydration shell of the Gdm+ central carbon,

236

which comprises about n ¯ C = 20.5 ± 0.5 water molecules. That value equals the number

237

of water molecules in the first hydration shell of the methane molecule. 39 The second peak

238

of gC (r) (located at about 6.3 ˚ A) is particularly flat. Lastly, gC (r) may be considered as

239

converged to ρ˜water as soon as 7 ˚ A. From the plots of the functions gC (r, 30◦ ) and gC∗ (r, 30◦ ),

240

it appears that the water molecules can be closer by about 0.2 ˚ A to the Gdm+ carbon atom

241

when located in the apical direction of the Gdm+ plane than when they are in the plane.

242

The function gN (r) is characterized by two sharp peaks located around 2.9 and 4.9 ˚ A.

11



243

However, their heights are only 1.60 and 1.25, respectively. Between these peaks, we note

244

a much flatter peak at about 3.7 ˚ A (its height is about 1.10). By comparing the profiles of

245

A of the function gN (r) the functions gN (r, 60◦ ) and gN∗ (r, 60◦ ), the two peaks at 2.9 and 4.9 ˚

246

correspond to water molecules interacting with the three Gdm+ NH2 moieties mainly in the

247

cation symmetry plane, while the gN (r) peak at 3.7 ˚ A arises mainly from water molecules

248

located above this plane. The weak heights of the gN (r) first peaks agree with all the

249

former experimental and theoretical conclusions regarding Gdm+ , i.e. the weakness of its

250

hydration structure, which is usually inferred to be at the origin of its protein denaturing

251

agent nature. 10

252

The main features of our radial distribution functions are in line with those computed

253

from the standard force field simulations (see Supporting Information). The main difference

254

between these two kinds of functions is the location of their first peaks. The standard force

255

field ones are shifted to larger distances compared to our polarizable model, by about 0.2

256

(nitrogen) and 0.5 (carbon) ˚ A for gX (r) and gX (r, φ0 ) functions. That was expected from

257

the comparison of quantum and standard force field structures of the dimer Gdm+ /H2 O (see

258

∗ Section 2.1). Regarding the complementary functions gX,φ , the polarizable and standard 0

259

force field ones remarkably agree, even if the first peak location of standard force field

260

functions is also shifted to larger Gdm+ /water distances.

261

In line with Warren and Patel, 40 we assume the sum of the water and Gdm+ molecular

262

radii Rwater +RGdm+ to correspond to the most probable carbon/water oxygen distance in the

263

Gdm+ first hydration shell. Our data yields then Rwater + RGdm+ = 3.8 ˚ A. From pure liquid

264

water simulations, the TCPE/2013 water model yields Rwater = 1.4 ˚ A. Hence RGdm+ = 2.4

265

˚ A, a value which is in line with that reported by Marcus, 2.1 ± 0.2 ˚ A. 9 The standard force

266

field leads to a slightly larger value for RGdm+ , about 2.6 ˚ A. Lastly, the data of our former

267

study focusing on methylated ammonium ions yields RNH+4 = 1.4 ˚ A.

12


Page 12 of 37

Page 13 of 37 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


268

3.2

269

3.2.1

270

The Gdm+ PMF profiles at droplet or at liquid water-vacuum interfaces and computed

271

from our polarizable approach are plotted as a function of the distance r or z in Figure 4,

272

+ together with the profiles of NH+ 4 and of K (only in the liquid water/vacuum case for the

273

latter cation). For the latter two ions their liquid water-vacuum profiles are taken from our

274

earlier study dealing with methylated ammonium cations. 25 The Gdm+ droplet PMF profiles

275

computed from the standard force field are provided as Supporting Information, and in the

276

liquid water/vacuum case, the PMF is also plotted in Figure 4(b) (together with the NH+ 4

277

one as computed also from a standard force field 36 ). The PMF profiles are set to zero at the

278

droplet centers and at the water slab cores.

279

280

Guanidinium properties at the water-vacuum interface Guanidinium propensity for water-vacuum interfaces

From the PMF functions, we define in the present study the ion propensity for the water/vacuum interface, PI, in the droplet case as ! Z r+ PMF(r) 1 Z r+ exp − P (r)dr = × 4πr2 dr, PI = Z r− kB T r−

(6)

281

here, P (r) is the ion location probability density, Z is the normalization factor, and r− and

282

r+ are the distances for which the water density drops from 95% to 5% of its value within

283

the droplet core (typically r+ − r− = 4 ˚ A). PI is thus the probability for an ion to be at the

284

droplet surface. In the liquid water/vacuum case, we define PI as

PI =

Z

z+

z−

h

!

i PMF(z) exp − − 1 dz, kB T

(7)

285

here, z − and z + are defined as above. In the slab case, the PI definition is reminiscent of

286

the function usually considered to estimate the surface excess/deficiency of a species. 41 The

287

PI numerical values are summarized in Figure 4.

288

In the droplet case, the PI value decreases as the droplet size increases regardless of

13



289

the ion and the force field. Our polarizable approach predicts the Gdm+ PI values to be

290

+ from 3 to 5 times smaller than the NH+ 4 ones, and the standard force field Gdm values are

291

in between the latter ones. In line with the standard force field, our polarizable approach

292

predicts thus Gdm+ to have a weaker propensity for the droplet surface than NH+ 4 . The plots

293

of Figure 4(b) also show that Gdm+ is predicted by both our our polarizable approach and

294

the standard force field to have a weaker propensity for the liquid water/vacuum interface

295

+ than both NH+ 4 and K .

296

Regarding our polarizable approach, the weakest propensity of Gdm+ for aqueous inter-

297

+ faces compared to NH+ 4 and K may be explained not only from both its large size and its

298

structure (that are inferred to be responsible for its hydration asymmetric properties 12,19 ),

299

but also because ion/water dispersion plays a stronger role for the Gdm+ hydration than for

300

+ the hydration of the two smaller cations NH+ 4 and K . From the plots of the mean ion/water

301

interaction energy values within the droplets (see Supporting Information), it appears that

302

the Gdm+ /water U¯ iw profile is very flat within the droplets, and that the ion/water energy

303

iw iw components U¯rep and U¯pol are centrifugal whereas the energy term U¯ ww and the ion/water

304

iw iw are centripetal. All these results are in line with our earlier and U¯qq energy components U¯disp

305

ones regarding both methylated ammonium and alkylated carboxylate ions. 25,28 They sug-

306

gest thus the hydration energetics of monovalent ions to obey a common pattern (within the

307

framework of our polarizable force field). However, the magnitude of the latter energy terms

308

iw values and components depends on the ion nature. In particular, the difference in the U¯disp

309

between the droplet core and the droplet surface is about twice as large for Gdm+ than for

310

25 + NH+ 4 (see Ref. ). Moreover, in the Gdm case, the latter difference is larger in magnitude

311

than the one corresponding to the sum of the electrostatic and polarization energy com-

312

ponents. As dispersion is an isotropic energy term favoring hydration within aqueous core

313

environment, that explains (at least partially) the weaker propensity of Gdm+ compared to

314

NH+ 4 for the droplet surface.

315

Another striking feature of the Gdm+ droplet PMF profiles computed from our polariz-

14


Page 14 of 37

Page 15 of 37 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


316

able approach is their very close properties among them and also as compared to the liquid

317

water-vacuum profile. That suggests the interfacial behavior of Gdm+ to be independent

318

from aqueous environments comprising at least 200 water molecules. That seems to be in

319

line with a recent infrared dissociation (IRPD) spectroscopy study of Heiles et al 15 of Gdm+

320

hydrated in small water environments (comprising from 5 to 100 water molecules) and show-

321

ing the profile of the IPRD spectra corresponding to NH bonded groups to be converged as

322

soon as Gdm+ hydrated aggregates comprising 50 water molecules.

323

However the weakest surface propensity of Gdm+ compared to NH+ 4 predicted by both

324

our polarizable approach and the standard force field contrasts with recent simulation find-

325

ings of Werner et all 19 agreeing with PhotoElectron spectroscopy experiments, i.e. Gdm+

326

− has a higher relative surface propensity than NH+ 4 in pure cation/Cl 2 M salt solutions.

327

Interestingly, the standard force field that we considered in the present study is the one used

328

by the latter authors. However we discuss in the present study the properties at the liquid

329

water-vacuum interface of hydrated cations in absence of explicit counter ions and at a much

330

smaller concentration (0.03 M) than in the latter studies. We may note also that simula-

331

tions based on the standard force field showed a dependence of the Gdm+ relative surface

332

propensity on the GdmCl salt concentration (in particular the Gdm+ surface propensity is

333

shown to increase with the salt concentration 18,19 ). Our results considered together with the

334

earlier ones seem thus to show the aqueous interfacial behavior of Gdm+ to depend strongly

335

on its concentration. However, as mentioned above, we don’t explicitly account for counter

336

ions in our interface simulations and halide ions like Cl− and Br− are considered to have a

337

strong propensity for water/vacuum interfaces. 42 Hence one may reasonably propose that

338

these anions can ”drive” cations to the surface. We are presently studying the interactions

339

of the latter halide ions with ammonium and guanidinium in water using our polarizable ap-

340

proach to investigate the role played by counter ions on the modulation of the guanidinium

341

propensity for the liquid water/vacuum interface.

15



Page 16 of 37

342

3.2.2

Guanidinium orientational properties at water-vacuum interfaces

343

Recent theoretical studies 17,18 concluded that a parallel orientation relative to the interface

344

of the Gdm+ symmetry plane is strongly favored close to the liquid water surface. Following

345

the ideas of Ou et al, 18 we investigate the Gdm+ orientation behavior by computing the mean

346

value of the second order Legendre polynomial along our droplet and liquid water-vacuum

347

simulations L2 (θ) =

1 3 cos2 θ − 1 , 2

(8)

348

where θ is the angle between the vector l⊥ and the vector connecting the droplet COM to the

349

Gdm+ carbon atom or between l⊥ and the axe orthogonal the liquid water-vacuum interface.

350

¯ 2 (dc ) computed along the umbrella sampling simulations The corresponding mean values L

351

based on our polarizable model are plotted as a function of the target distances dc in Figure

352

¯ 2 ≈ 0) from the droplet or 5. From that plot, the Gdm+ plane orientation is random (i.e. L

353

the water slab core up to about 3 ˚ A before the interface, whereas a parallel orientation is

354

¯ 2 ≥ 0.7), regardless of the water system (droplet or bulk). favored close to the interface (L

355

The transition in the Gdm+ orientation is particularly abrupt (it is achieved 1 ˚ A before the

356

interface).

357

Regarding the liquid water-vacuum interface, our result agrees with simulations generated

358

using either a standard pairwise force field 17 and a fluctuating charge polarizable approach. 18

359

Hence a parallel orientation of Gdm+ to the liquid water-vacuum interface is predicted by

360

all the MD simulations reported to date, regardless of the force field sophistication and the

361

parameterization protocol used to assign their parameters. That result clearly contrasts with

362

a recent heterodyne-detected vibrational sum frequency generation study of Strazdaite et al 43

363

who concluded a near-parallel orientation to the liquid water-air surface to be observed for

364

only a minor fraction of the methyl guanidinium ions.

16


Page 17 of 37 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


365

3.3

Guanidinium bulk hydration enthalpy

366

As in our recent studies, 25,26,28 we extrapolated single ion bulk hydration enthalpies by

367

adjusting the parameters of the following power law function of the droplet size N

iw f (N ) = U∞ +

a N 1/3

(9)

368

to reproduce the mean ion/water potential energy in droplets U¯ iw (N ) computed from the

369

umbrella sampling simulations according to

U¯ iw (N ) =

X

U¯ iw (rci ) × P (rci ).

(10)

i

370

Here the sum runs on the umbrella sampling simulations whose target distances are rci .

371

ww Denoting by ∆Ubulk the bulk solvent destabilization energy induced by the ion presence, we

372

iw ww assumed in our former studies that the sum U∞ +∆Ubulk yields the ion hydration microscopic

373

iw + energy in bulk water and thus the ion bulk hydration enthalpy according to ∆Hhyd = U∞

374

ww − RT . Below, if not discussed, the error estimates corresponding to mean energy ∆Ubulk

375

values are computed from the root mean square deviations of the data extracted along

376

the last 9 ns segment of the MD simulations and by considering these data as temporally

377

uncorrelated (along each simulations, 9 000 regularly spaced data are extracted).

378

As several authors, 15,44–46 we thus extrapolate bulk data from isolated finite size system

379

results. A priori, the reliability of such an extrapolation process is physically founded if the

380

finite size systems are simulated by accounting for bulk environment effects on them. 47,48

381

However we may note that from our own TCPE/2013 computations the interfacial potential

382

at the water-vacuum boundary in droplets whose size is ≥ 100 is already converged to the

383

bulk value (this will be discussed elsewhere). That suggests that the water HB network struc-

384

ture at the water-vacuum interface in intermediate size droplets is predicted by the model

385

TCPE/2013 to be close to the bulk one and, thus, that justifies a priori our extrapolation

17



386

protocol.

387

ww is computed by comparing the total water interaction energy from a Gdm+ 10 ∆Ubulk

388

ns bulk simulation to that of a pure water one (and comprising the same number of water

389

ww molecules as for Gdm+ ). For our polarizable model, ∆Ubulk = +29.7 ± 0.8 kcal mol−1 , a

390

25 value which is smaller by about 6 kcal mol−1 than for NH+ Hence, even if Gdm+ is a 4.

391

˚ larger ion than NH+ 4 (their molecular radii are 2.4 and 1.4 A, respectively), it destabilizes

392

less the bulk water HB network than NH+ 4 . That can be a consequence of the more diluted

393

+ ion charge in Gdm+ than in NH+ 4 and of the absence of hydrophobic faces/moieties in NH4 .

394

ww value for Gdm+ is 40.1 ± 1.1 kcal Computed accordingly, the standard force field ∆Ubulk

395

mol−1 , a value larger by 25 % than the polarizable model estimate.

396

In Figure 6, we plot the values U¯ iw (N ) as a function of N −1/3 for both Gdm+ and NH+ 4

397

solvated in the four water droplets considered in the present study, as well as the results of

398

the linear regression process to adjust the f (N ) parameters (in that figure, we also reported

399

iw the data for Gdm+ computed from the standard force field). For Gdm+ , the parameter U∞

400

computed from our polarizable approach is -112.6 kcal mol−1 , a value 10 kcal mol−1 smaller

401

+ −1 iw than for NH+ 4 (-122.4 kcal mol ). Note that if we extrapolate the value U∞ for NH4 by

402

considering all together our earlier droplet data 25,26 and the present new ones (that leads to

403

a data set including the results of droplets whose size ranges from 50 up to 10 000 molecules),

404

we get then again the same value than the latter one, within 0.1 kcal mol−1 .

405

As discussed in Section 3.2.1, the U¯ iw (rci ) profiles are very flat, regardless of the ion loca-

406

tion within the droplet. As already discussed for both methylated ammonium and alkylated

407

carboxylate ions, 26 the values U¯ iw (N) depend thus marginally on the features of the ion

408

location probability density P (r). As for the latter two kinds of ions, we estimate here also

409

that the uncertainty affecting the values U¯ iw (N ) and arising from possible P (r) drawbacks

410

iw to be 1 kcal mol−1 . As the uncertainty regarding the parameters U∞ and originating from

411

the linear fitting process is 0.3 kcal mol−1 and the uncertainty affecting the values U¯ iw (N)

412

arising from the finite size of the data extracted along the MD simulations is 0.1 kcal mol−1 ,

18


Page 18 of 37

Page 19 of 37 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

413


iw the total uncertainty affecting the parameters U∞ is 1.4 kcal mol−1 .

414

iw ww computed from our polarizable approach and U∞ The above results regarding ∆Ubulk

415

are in line with the common assumption that the Gdm+ /water interactions are weaker than

416

+ the NH+ 4 ones. They lead to a Gdm bulk hydration enthalpy ∆Hhyd of -82.9 ± 2.2 kcal

417

mol−1 , which is about 4 kcal mol−1 weaker than the NH+ 4 one computed using the same

418

theoretical protocol (-87.1 ± 2.2 kcal mol−1 ). Our Gdm+ estimate is almost half as large

419

as the experimental-based one reported recently by Marcus, 9 about -143.8 ± 1.9 kcal mol−1

420

(however, that value is considered as puzzling by the latter author). Assuming the reliability

421

of the hydration entropy component −T∆Shyd computed from the data reported by Mar-

422

cus 9 for Gdm+ , 4.5±0.5 kcal mol−1 at 300 K, that yields 78.4±2.6 kcal mol−1 for the single

423

hydration free energy of Gdm+ , ∆Ghyd (Gdm+ ). Regarding NH+ 4 , by considering the hydra-

424

tion entropy component −T∆Shyd value of 6.6 kcal mol−1 reported in Ref., 14 our theoretical

425

−1 protocol yields thus ∆Ghyd (NH+ 4 )= −80.5 ± 2.2 kcal mol , in very good agreement with

426

the experimental-based estimate, about -81 kcal mol−1 . 14 We computed also the absolute

427

bulk hydration free energies ∆Ghyd for Gdm+ and NH+ 4 using our polarizable approach and a

428

standard Thermodynamical Integration scheme (by considering 40 intermediate steps and by

429

accounting for all the corrections arising from periodic condition drawbacks, like the absence

430

of explicit water/vacuum interface). For the two ions, we get then ∆Ghyd values very close

431

to that above ones, i.e. -76 and -78 kcal mol−1 for Gdm+ and NH+ 4 , respectively, however

432

with a large uncertainty for Gdm+ , about ±5 kcal mol−1 a value three times larger than

433

+ for NH+ 4 (this arises from the very large intra molecular electrostatic energy within Gdm

434

for our polarizable approach, about 250 kcal mol−1 , and because our TI decoupling scheme

435

scales all the atomic static charges and not only the Gdm+ /water electrostatic interactions).

436

iw Computed accordingly, the standard force field U∞ value is -130.6 ± 1.4 kcal mol−1 ,

437

larger by 18 kcal mol−1 than the one computed by our polarizable approach. Accounting for

438

ww ∆Ubulk , the standard force field yields a ∆Hhyd value of -90.5 ± 1.9 kcal mol−1 and a ∆Ghyd

439

value of -86.0 kcal mol−1 for Gdm+ . Both the latter value are 10 kcal mol−1 larger than our

19



440

polarizable model estimates. However the standard force field ∆Hhyd (Gdm+ ) value is also

441

largely weaker than the experimental estimate reported by Marcus, by about 54 kcal mol−1 .

442

The ∆Ghyd (Gdm+ ) value computed by our polarizable approach is by 18 kcal mol−1

443

larger in magnitude than the ones computed from quantum computations coupled with a

444

polarizable continuous solvent approach (QM/PCM), 13 about 58 kcal mol−1 . However, the

445

results of this kind of coupled theoretical approach are highly sensitive to a set of empir-

446

ical parameters, the atomic radii, whose quality is assessed by comparing to experiment

447

the QM/PCM estimates regarding the experimental bulk hydration free energies of a set of

448

training molecules and ions, 49,50 and as far as we know, no reliable experimental thermody-

449

namical data is available for Gdm+ . It is thus not obvious to comment here the difference

450

between our ∆Ghyd (Gdm+ ) estimate and the quantum ones. Nevertheless, we may note

451

that the QM/PCM ∆Ghyd value for Gdm+ is also largely underestimated compared to that

452

of the methyl-guanidinium ion and computed using the popular macroscopic solvent model

453

FDPB/γ developed by Honig and co-workers, 51 -66.1 kcal mol−1 .

454

Considering our polarizable approach results, the bulk hydration enthalpy and free energy

455

values for Gdm+ are only a few kcal mol−1 smaller than the NH+ 4 ones. This arises from

456

ww iw discussed above and ∆Ubulk two phenomena with opposite energetic trends as the values U∞

457

suggest. First, the Gdm+ /water interactions are weaker than NH+ 4 ones, which in turn leads

458

to a weaker perturbation of the bulk water HB network in presence of Gdm+ as compared

459

+ to NH+ 4 . In all, this is in line with the well accepted property of Gdm , i.e. it interacts

460

“weakly” with water.

461

4

462

We report a theoretical analysis of the hydration process of the Gdm+ in finite and infinite

463

aqueous environments from microscopic simulations performed using a polarizable force field

464

approach. Our results show Gdm+ to have a strong propensity for the core of aqueous

Conclusions

20


Page 20 of 37

Page 21 of 37 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


465

environment and a weak propensity for water/vacuum interfaces. In particular, the Gdm+

466

propensity for water-vacuum interface is predicted by our approach to be weaker than the

467

+ NH+ 4 and K ones. This a priori conflicts with earlier simulations based on standard pairwise

468

force field approaches. 17,18 However, besides theoretical methodology differences, it has to

469

be noted here that we simulate the hydration process of Gdm+ in bulk water at an ion

470

concentration that is about two order of magnitude smaller than in the latter theoretical

471

studies (0.03 and 2 to 5 M, respectively). Nevertheless, as already theoretically shown, 17,18

472

our approach predicts a preferential parallel orientation of the Gdm+ symmetry plane to

473

the water surface near water-vacuum interfaces, which is in contradiction with a recent

474

experimental study. 43

475

For both the cations Gdm+ and NH+ 4 , we estimate their absolute bulk hydration en-

476

thalpy ∆Hhyd and free energy ∆Ghyd values by extrapolating droplet data computed from

477

our polarizable model and by considering available hydration entropy values. Our NH+ 4

478

thermodynamic values match the experimental-based ones, within the uncertainty bars. In

479

−1 + particular, we get ∆Ghyd (NH+ 4 ) = -80.1 ± 2.2 kcal mol . Regarding Gdm , both our ∆Hhyd

480

and ∆Ghyd values are a few kcal mol−1 smaller in magnitude than their NH+ 4 counterparts,

481

i.e. ∆Hhyd (Gdm+ ) = -82.9 ± 2.2 kcal mol−1 and ∆Ghyd (Gdm+ )= -78.4 ± 2.6 kcal mol−1 .

482

Compared to monoatomic alkali cations, the latter difference between the Gdm+ and NH+ 4

483

∆Ghyd estimates is close to that existing between Rb+ and K+ . 52

484

Our ∆Hhyd (Gdm+ ) value is about half as large as the experimental based one reported by

485

Marcus 9 that is clearly puzzling as it is even larger in magnitude than all the available ∆Hhyd

486

estimates for the smallest alkali cation Li+ by at least 10 kcal mol−1 . On the other hand,

487

our ∆Ghyd (Gdm+ ) value is by about 18 kcal mol−1 larger in magnitude than theoretical

488

quantum based estimates, 13 which can be considered as surprisingly low when compared

489

to the target hydration free energy value of the methyl guanidinium cation considered in

490

popular theoretical methods like the FDPB/γ one developed by Honig and co-workers 51

491

(-58 and -66 kcal mol−1 , respectively). Within our polarizable theoretical framework, the

21



492

overall large hydration energy value of Gdm+ arises mainly from the weak perturbation of

493

the water bulk HB network due to the Gdm+ presence, that compensates weak Gdm+ /water

494

interactions (at least weaker than for NH+ 4 ).

495

For comparison purposes, we also computed the Gdm+ single ion hydration thermody-

496

namic data from the same theoretical protocol, however, by considering the standard (i.e.

497

non polarizable) force field that was shown by Wernesson et al to provide accurate results

498

regarding the Gdm+ salts behavior at the liquid water interface. If most of our polarizable

499

model predictions agree qualitatively with the standard force field ones, we may note that

500

the standard force field leads to overestimate both the Gdm+ single ion hydration enthalpy

501

and free energy values compared to our polarizable approach, by about 10 kcal mol−1 . In

502

all, a standard force field approach may be thus considered as providing an accurate enough

503

picture of the hydration process of complex organic ions like Gdm+ .

504

The relatively large magnitude of our hydration energy values for Gdm+ suggests a

505

relatively high affinity of that cation for aqueous environment. Hence the Gdm+ original

506

protein denaturing properties have to result from a competition between Gdm+ hydration

507

effects and Gdm+ /solute interactions that can be modulated by weak differences in the

508

protein chemical environment or by the nature of the counter ion with which it is paired. 53

509

That suggests in particular the pivotal role played by the Gdm+ hydrophobic faces to explain

510

the original properties of that cation rather than a weak propensity for aqueous environments.

511

Supporting Information Available

512

This material is available free of charge via the Internet at http://pubs.acs.org. It pro-

513

vides more detailed discussions about the polarizable force field accuracy to model gas

514

phase hydrated Gdm+ clusters, the force field parameters, the quantum VDZ/VTZ/VQZ

515

Gdm+ /water binding energies and a comparison with the data of both the standard and

516

the polarizable force field, the radial elongation energy profiles of a single water molecule

22


Page 22 of 37

Page 23 of 37 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


517

interacting with Gdm+ within its symmetry plane or orthogonally to it, the stepwise binding

518

energies for the Gdm+ micro hydrated structures where a water molecule is in the cation sec-

519

ond hydration shell, the angular distribution function corresponding to the angle CN − Owater

520

in bulk water, a comparison of the radial Gdm+ /water distribution functions computed from

521

standard and polarizable force field bulk simulations, the Gdm+ location probability densi-

522

ties within the droplets and the Gdm+ /water interaction energy decomposition in droplets

523

(from the polarizable model).

524

Acknowledgments

525

This work was granted access to the HPC resources of [CCRT/CINES/IDRIS] under the

526

allocation 2016-[6100] by GENCI (Grand Equipement National de Calcul Intensif).

527

References

528

(1) Greene, R. F.; Pace, C. N. Urea and Guanidine Hydrochloride Denaturation of Ribonu-

529

clease, Lysozyme, Alpha-Chymotrypsin, and Beta-Lactoglobulin. Journal of Biological

530

Chemistry 1974, 249, 5388–5393.

531

(2) Okur, H. I.; Hladlkov, J.; Rembert, K. B.; Cho, Y.; Heyda, J.; Dzubiella, J.; Cre-

532

mer, P. S.; Jungwirth, P. Beyond the Hofmeister Series: Ion-Specific Effects on Pro-

533

teins and Their Biological Functions. The Journal of Physical Chemistry B 2017, 121,

534

1997–2014.

535

(3) Shih, O.; England, A. H.; Dallinger, G. C.; Smith, J. W.; Duffey, K. C.; Cohen, R. C.;

536

Prendergast, D.; Saykally, R. J. Cation-Cation Contact Pairing in Water: Guanidinium.

537

The Journal of Chemical Physics 2013, 139, 035104.

538

(4) Mason, P. E.; Neilson, G. W.; Enderby, J. E.; Saboungi, M.-L.; Dempsey, C. E.; MacK-

23



539

erell, A. D.; Brady, J. W. The Structure of Aqueous Guanidinium Chloride Solutions.

540

Journal of the American Chemical Society 2004, 126, 11462–11470.

541

(5) Kub´ıˇcková, A.; Kˇr´ıˇzek, T.; Coufal, P.; Wernersson, E.; Heyda, J.; Jungwirth, P. Guani-

542

dinium Cations Pair with Positively Charged Arginine Side Chains in Water. The Jour-

543

nal of Physical Chemistry Letters 2011, 2, 1387–1389.

544

(6) Magalhaes, A.; Maigret, B.; Hoflack, J.; Gomes, N. F.; Scheraga, H. Contribution of

545

Unusual Arginine-Arginine Short-Range Interactions to Stabilization and Recognition

546

in Proteins. Journal of Protein Chemistry 1994, 13, 195–215.

547

(7) Vondráˇsek, J.; Mason, P. E.; Heyda, J.; Collins, K. D.; Jungwirth, P. The Molecular

548

Origin of Like-Charge Arginine-Arginine Pairing in Water. The Journal of Physical

549

Chemistry B 2009, 113, 9041–9045.

550

(8) Mason, P. E.; Dempsey, C. E.; Neilson, G. W.; Kline, S. R.; Brady, J. W. Preferential

551

Interactions of Guanidinum Ions with Aromatic Groups over Aliphatic Groups. Journal

552

of the American Chemical Society 2009, 131, 16689–16696.

553

554

(9) Marcus, Y. The Guanidinium Ion. The Journal of Chemical Thermodynamics 2012, 48, 70 – 74.

555

(10) Mason, P. E.; Neilson, G. W.; Dempsey, C. E.; Barnes, A. C.; Cruickshank, J. M. The

556

Hydration Structure of Guanidinium and Thiocyanate Ions: Implications for Protein

557

Stability in Aqueous Solution. Proceedings of the National Academy of Sciences 2003,

558

100, 4557–4561.

559

(11) Scott, J. N.; Nucci, N. V.; Vanderkooi, J. M. Changes in Water Structure Induced by

560

the Guanidinium Cation and Implications for Protein Denaturation. The Journal of

561

Physical Chemistry A 2008, 112, 10939–10948.

24


Page 24 of 37

Page 25 of 37 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


562

(12) Cooper, R. J.; Heiles, S.; DiTucci, M. J.; Williams, E. R. Hydration of Guanidinium:

563

Second Shell Formation at Small Cluster Size. The Journal of Physical Chemistry A

564

2014, 118, 5657–5666.

565

(13) Ciftcioglu, G. A.; Trindle, C. Computational Estimates of Thermochemistry and pKa

566

Values of Cyclopropenyl Imine Superbases. International Journal of Quantum Chem-

567

istry 2014, 114, 392–399.

568

(14) Floriiàn, J.; Warshel, A. Calculations of Hydration Entropies of Hydrophobic, Polar,

569

and Ionic Solutes in the Framework of the Langevin Dipoles Solvation Model. The

570

Journal of Physical Chemistry B 1999, 103, 10282–10288.

571

572

(15) Heiles, S.; Cooper, R. J.; DiTucci, M. J.; Williams, E. R. Hydration of Guanidinium Depends on its Local Environment. Chemical Science 2015, 6, 3420–3429.

573

(16) Mountain, R. D.; Thirumalai, D. Alterations in Water Structure Induced by Guani-

574

dinium and Sodium Ions. The Journal of Physical Chemistry B 2004, 108, 19711–

575

19716.

576

(17) Wernersson, E.; Heyda, J.; Vazdar, M.; Lund, M.; Mason, P. E.; Jungwirth, P. Ori-

577

entational Dependence of the Affinity of Guanidinium Ions to the Water Surface. The

578


579

(18) Ou, S.; Lucas, T. R.; Zhong, Y.; Bauer, B. A.; Hu, Y.; Patel, S. Free Energetics and

580

the Role of Water in the Permeation of Methyl Guanidinium across the BilayerWater

581

Interface: Insights from Molecular Dynamics Simulations Using Charge Equilibration

582

Potentials. The Journal of Physical Chemistry B 2013, 117, 3578–3592.

583

¨ (19) Werner, J.; Wernersson, E.; Ekholm, V.; Ottosson, N.; Ohrwall, G.; Heyda, J.; Pers-

584

son, I.; Söderström, J.; Jungwirth, P.; Björneholm, O. Surface Behavior of Hydrated

585

Guanidinium and Ammonium Ions: A Comparative Study by Photoelectron Spec-

25



586

troscopy and Molecular Dynamics. The Journal of Physical Chemistry B 2014, 118,

587

7119–7127.

588

(20) Vazdar, M.; Uhlig, F.; Jungwirth, P. Like-Charge Ion Pairing in Water: An Ab Initio

589

Molecular Dynamics Study of Aqueous Guanidinium Cations. The Journal of Physical

590

Chemistry Letters 2012, 3, 2021–2024.

591

(21) Vorobyov, I.; Li, L.; Allen, T. W. Assessing Atomistic and Coarse-Grained Force Fields

592

for Protein Lipid Interactions: the Formidable Challenge of an Ionizable Side Chain in

593

a Membrane. The Journal of Physical Chemistry B 2008, 112, 9588–9602.

594

(22) Wu, J. C.; Chattree, G.; Ren, P. Automation of AMOEBA Polarizable Force Field

595

Parameterization for Small Molecules. Theoretical Chemistry Accounts 2012, 131, 1138.

596

(23) Shi, Y.; Xia, Z.; Zhang, J.; Best, R.; Wu, C.; Ponder, J. W.; Ren, P. Polarizable Atomic

597

Multipole-Based AMOEBA Force Field for Proteins. Journal of Chemical Theory and

598

Computation 2013, 9, 4046–4063.

599

(24) Réal, F.; Vallet, V.; Flament, J.-P.; Masella, M. Revisiting a Many-Body Model for

600

Water Based on a Single Polarizable Site. From Gas Phase Clusters to Liquid and

601

Air/Liquid Water Systems. The Journal of Chemical Physics 2013, 139, 114502.

602

(25) Houriez, C.; Meot-Ner (Mautner), M.; Masella, M. Simulated Solvation of Organic Ions:

603

Protonated Methylamines in Water Nanodroplets. Convergence toward Bulk Properties

604

and the Absolute Proton Solvation Enthalpy. The Journal of Physical Chemistry B

605

2014, 118, 6222–6233.

606

(26) Coles, J. P.; Houriez, C.; Meot-Ner (Mautner), M.; Masella, M. Extrapolating Sin-

607

gle Organic Ion Solvation Thermochemistry from Simulated Water Nanodroplets. The

608


26


Page 26 of 37

Page 27 of 37 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


609

(27) Réal, F.; Gomes, A. S. P.; Mart´ınez, Y. O. G.; Ayed, T.; Galland, N.; Masella, M.;

610

Vallet, V. Structural, Dynamical, and Transport Properties of the Hydrated Halides:

611

How Do At− bulk Properties Compare with those of the Other Halides, from F− to I− ?

612

The Journal of Chemical Physics 2016, 144, 124513.

613

(28) Houriez, C.; Meot-Ner (Mautner), M.; Masella, M. Simulated Solvation of Organic Ions

614

II: Study of Linear Alkylated Carboxylate Ions in Water Nanodrops and in Liquid Wa-

615

ter. Propensity for Air/Water Interface and Convergence to Bulk Solvation Properties.

616

The Journal of Physical Chemistry B 2015, 119, 12094–12107.

617

(29) Meot-Ner, M. Heats of Hydration of Organic ions: Predictive Relations and Analysis of

618

Solvation Factors Based on Ion Clustering. The Journal of Physical Chemistry 1987,

619

91, 417–426.

620

621

(30) Meot-Ner (Mautner), M. The Ionic Hydrogen Bond. Chemical Reviews 2005, 105, 213– 284.

622

(31) Tuttle, T. R.; Malaxos, S.; Coe, J. V. A New Cluster Pair Method of Determining Ab-

623

solute Single Ion Solvation Energies Demonstrated in Water and Applied to Ammonia.

624

The Journal of Physical Chemistry A 2002, 106, 925–932.

625

(32) Frisch, M.; Trucks, G.; Schlegel, H.; Scuseria, G.; Robb, M. A.; Cheeseman, J. R.;

626

Scalmani, G.; Barone, V.; Mennucci, B.; Petersson, G. A. et al. Gaussian 09 Revision

627

D.01. Gaussian Inc. Wallingford CT 2009.

628

(33) Feller, D. The Use of Systematic Sequences of Wave Functions for Estimating the

629

Complete Basis Set, Full Configuration Interaction Limit in Water. The Journal of

630

Chemical Physics 1993, 98, 7059–7071.

631

(34) Helgaker, T.; Klopper, W.; Koch, H.; Noga, J. Basis-Set Convergence of Correlated

632

Calculations on Water. The Journal of Chemical Physics 1997, 106, 9639–9646.

27



633

(35) Brooks, B.; Brooks III, C.; Mackerell, A.; Nilsson, L.; Petrella, R.; Roux, B.; Won, Y.;

634

Archontis, G.; Bartels, C.; Boresch, S. et al. CHARMM: The Biomolecular Simulation

635

Program. Journal of Computational Chemistry 2009, 30, 1545–1615.

636

637

(36) Boudon, S.; Wipff, G. Free Energy Calculations Involving NH4+ in Water. Journal of Computational Chemistry 1991, 12, 42–51.

638

(37) Toukmaji, A.; Sagui, C.; Borad, J.; Darden, T. Efficient Particle-Mesh Ewald Based

639

Approach to Fixed and Induced Dipolar Interactions. The Journal of Chemical Physics

640

2000, 113, 10913–10927.

641

(38) Kästner, J.; Thiel, W. Bridging the Gap between Thermodynamic Integration and Um-

642

brella Sampling Provides a Novel Analysis Method: Umbrella Integration. The Journal

643

of Chemical Physics 2005, 123, 144104.

644

(39) Asthagiri, D.; Merchant, S.; Pratt, L. R. Role of Attractive Methane-Water Interactions

645

in the Potential of Mean Force between Methane Molecules in Water. The Journal of

646

Chemical Physics 2008, 128, 244512.

647

(40) Warren, G. L.; Patel, S. Hydration Free Energies of Monovalent Ions in Transferable

648

Intermolecular Potential Four Point Fluctuating Charge Water: An Assessment of Sim-

649

ulation Methodology and Force Field Performance and Transferability. The Journal of

650

Chemical Physics 2007, 127, 064509.

651

652

653

654

655

(41) Manciu, M.; Ruckenstein, E. Specific Ion Effects via Ion Hydration: I. Surface Tension. Advances in Colloid and Interface Science 2003, 105, 63 – 101. (42) Piatkowski, L.; Zhang, Z.; Backus, E. H. G.; Bakker, H. J.; Bonn, M. Extreme Surface Propensity of Halide Ions in Water. Nature Communications 2014, 5, 4083. (43) Strazdaite, S.; Versluis, J.; Ottosson, N.; Bakker, H. J. Orientation of Methylguani-

28


Page 28 of 37

Page 29 of 37 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


656

dinium Ions at the Water-Air Interface. The Journal of Physical Chemistry C 0, 0,

657

null.

658

(44) Donald, W. A.; Leib, R. D.; Demireva, M.; O’Brien, J. T.; Prell, J. S.; Williams, E. R.

659

Directly Relating Reduction Energies of Gaseous Eu(H2 O)3+ n , n = 55-140, to Aqueous

660

Solution: The Absolute SHE Potential and Real Proton Solvation Energy. Journal of

661

the American Chemical Society 2009, 131, 13328–13337.

662

(45) Coe, J. V. Connecting Cluster Anion Properties to Bulk: Ion Solvation Free Energy

663

Trends with Cluster Size and the Surface vs Internal Nature of Iodide in Water Clusters.

664

The Journal of Physical Chemistry A 1997, 101, 2055–2063.

665

(46) Coe, J. V.; Williams, S. M.; Bowen, K. H. Photoelectron Spectra of Hydrated Elec-

666

tron Clusters vs. Cluster Size: Connecting to Bulk. International Reviews in Physical

667

Chemistry 2008, 27, 27–51.

668

(47) Beglov, D.; Roux, B. Finite Representation of an Infinite Bulk System: Solvent Bound-

669

ary Potential for Computer Simulations. The Journal of Chemical Physics 1994, 100,

670

9050–9063.

671

(48) Darden, T.; Pearlman, D.; Pedersen, L. G. Ionic Charging Free Energies: Spherical

672

Versus Periodic Boundary Conditions. The Journal of Chemical Physics 1998, 109,

673

10921–10935.

674

675

676

677

678

(49) Cramer, C. J.; Truhlar, D. G. Implicit Solvation Models: Equilibria, Structure, Spectra, and Dynamics. Chemical Reviews 1999, 99, 2161–2200. (50) Tomasi, J.; Mennucci, B.; Cammi, R. Quantum Mechanical Continuum Solvation Models. Chemical Reviews 2005, 105, 2999–3094. (51) Sitkoff, D.; Sharp, K. A.; Honig, B. Accurate Calculation of Hydration Free Energies

29



679

Using Macroscopic Solvent Models. The Journal of Physical Chemistry 1994, 98, 1978–

680

1988.

681

(52) Lamoureux, G.; Roux, B. Absolute Hydration Free Energy Scale for Alkali and Halide

682

Ions Established from Simulations with a Polarizable Force Field. The Journal of Phys-

683

ical Chemistry B 2006, 110, 3308–3322.

684

(53) Heyda, J.; Okur, H. I.; Hlad´ılková, J.; Rembert, K. B.; Hunn, W.; Yang, T.; Dzu-

685

biella, J.; Jungwirth, P.; Cremer, P. S. Guanidinium Can both Cause and Prevent the

686

Hydrophobic Collapse of Biomacromolecules. Journal of the American Chemical Society

687

2017, 139, 863–870.

30


Page 30 of 37

Page 31 of 37 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


Figure 1: Most stable Gdm+ /(H2 O)n aggregate structures for n = 1 − 3 from MP2 quantum computations, from our polarizable and a standard force field. 19 The energies provided correspond to the stepwise binding energies ∆Un−1,n . The quantum data are provided in brackets and the standard force field ones in parentheses.

31



Figure 2: A water molecule interacting orthogonally to the Gdm+ plane.

32


Page 32 of 37

Page 33 of 37 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


Figure 3: Radial distribution functions centered at the Gdm+ nitrogen (a) and carbon (b) atoms from our polarizable force field. Black line: full functions, green lines : functions gX (r, φ0 ), blue line : complementary functions gX∗ (r, φ0 ). The nitrogen functions are averaged over the data of the three Gdm+ nitrogen atoms.

33



Figure 4: Ion/water PMFs at the water-vacuum interfaces. (a) Droplet profiles, in bold and dashed lines, PMF profiles for Gdm+ and NH+ 4 , respectively, the results for 200, 400, 600 and 800 water droplets are shown in black, blue, green and orange, respectively. (b) + 25 liquid water-vacuum profiles, black: Gdm+ , blue: NH+ 4 , red: K profile (from Ref. ) (bold 19 lines, our polarizable model, and dashed lines, the standard force field ). The vertical dashed lines show the water-vacuum boundaries. The numerical values correspond to the ion interface propensities PI. For droplets, and they are expressed in % (in parentheses, the ˚ NH+ 4 ones); and for the liquid water-vacuum interface, they are expressed in A (in parentheses the standard force field values).

34


Page 34 of 37

Page 35 of 37 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


Figure 5: Gdm+ orientation at water-vacuum interfaces. (a) droplet and (b) liquid watervacuum case from our polarizable force field. θ is the angle between the vector l⊥ and the vector connecting the droplet COM to the Gdm+ carbon atom (a) or between l⊥ and the vector orthogonal to the liquid water-vacuum interface (b). The results for 200, 400, 600 and 800 water molecule droplets are shown in black, blue, green and orange, respectively.

35



Figure 6: Gdm+ (black: our polarizable force-field approach, grey: standard force field −1/3 ¯ iw approach 19 ) and NH+ , 4 (blue) - water mean interaction energy U (N ) as a function of N where N is the water droplet size. Full symbols correspond to the droplets considered in the 25,26 present study. For NH+ regarding 4 , the empty symbols correspond to our earlier data droplets whose size varies from N = 50 up to 10 000. Dashed lines: results of the linear regression (the regression coefficients are greater than 0.99, regardless of the ions and the force field used).

36


Page 36 of 37

Page 37 of 37 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24



Solvation of the Guanidinium Ion in Pure Aqueous Environments: A

Recommend Documents