Interaction Energies in Complexes of Zn and Amino Acids: A

Subscriber access provided by UB + Fachbibliothek Chemie | (FU-Bibliothekssystem)

Article

Interaction Energies in Complexes of Zn and Amino-Acids: A Comparison of Ab-Initio and Force Field based Calculations Ran Friedman, Emma Ahlstrand, and Kersti Hermansson J. Phys. Chem. A, Just Accepted Manuscript • DOI: 10.1021/acs.jpca.6b12969 • Publication Date (Web): 08 Mar 2017 Downloaded from http://pubs.acs.org on March 14, 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 A 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 39

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

Interaction Energies in Complexes of Zn and Amino-Acids: a Comparison of Ab-Initio and Force Field based Calculations Emma Ahlstrand,†,‡ Kersti Hermansson,¶ and Ran Friedman∗,†,‡ 1

†Department of Chemistry and Biomedical Sciences, Linnæus University, 391 82 Kalmar, Sweden ‡Linnæus University Centre for Biomaterials Chemistry ¶Department of Chemistry, Ångström Laboratory, Uppsala University, BOX 538, 751 21 Uppsala, Sweden E-mail: [email protected] Phone: +46 (0)480 446290

1 ACS Paragon Plus Environment


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

2

3

Zinc plays important roles in structural stabilization of proteins, enzyme catalysis

4

and signal transduction. Many Zn binding sites are located at the interface between

5

the protein and the cellular fluid. In aqueous solutions, Zn ions adopt an octahedral

6

coordination, while in proteins zinc can have different coordinations with a tetrahedral

7

conformation found most frequently. The dynamics of Zn binding to proteins, and the

8

formation of complexes that involve Zn, are dictated by interactions between Zn and

9

its binding partners. We calculated the interaction energies between Zn and ligands in

10

complexes that mimic protein binding sites, and in Zn-complexes of water and one or

11

two amino acid moieties, using quantum-mechanics (QM) and molecular mechanics

12

(MM). It was found that MM calculations that neglect or only approximate polarizability

13

did not reproduce even the relative order of the interaction energies in these com-

14

plexes. Interaction energies calculated with the CHARMM-Drude polarizable force

15

field agreed better with the ab initio results, although the deviations between QM and

16

MM were still rather large (40–96 kcal/mol). In order to gain further insight into Zn-

17

ligand interactions, the free energies of interaction were estimated by QM calculations

18

with continuum solvent representation, and we performed energy decomposition anal-

19

ysis calculations to examine the characteristics of the different complexes. The ligands

20

were found to have high impact on the relative strength of polarization and electrostatic

21

interactions. Interestingly, ligand-ligand interactions did not play a significant role in

22

the binding of Zn. Finally, analysis of ligand exchange energies suggests that car-

23

boxylates could be exchanged with water molecules which explains the flexibility in Zn

24

binding dynamics. An exchange between carboxylate (Asp/Glu) and imidazole (His)

25

is less likely.

26

Introduction

27

Zn2+ is one of the most commonly encountered protein metal cofactors and the interplay

28

between Zn2+ and amino acid residues therefore plays an important role in biology. The 2 ACS Paragon Plus Environment

Page 2 of 39

Page 3 of 39

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


29

zinc ion participates in enzyme catalysis, structural stabilization of protein motifs, and bio-

30

logical regulation in enzymes. 1 Many of the catalytic Zn-binding sites are partially exposed

31

to solvent (intra- or extracellular fluid). The interactions between Zn and its ligands in the

32

protein ensure that the formation of Zn-protein complexes is thermodynamically and ki-

33

netically favored over its fully hydrated conformation. Moreover, the Zn-ligand interactions

34

tune the subtle dynamics and electrostatics of the protein. A sound description of the ion

35

in both these environments, i.e. protein and solution, is required for an understanding of

36

the ion’s binding and dynamics. In dilute aqueous solution, a Zn2+ ion is known to adopt

37

an octahedral coordination by binding to the oxygen atoms of six water molecules. 2,3 In

38

proteins, on the other hand, a tetrahedral conformation is most frequently found, 4 and the

39

most common amino acid ligands are cysteine, histidine, aspartate and glutamate. 5–8 The

40

ability of zinc to readily adopt different coordinations is an important feature of this ion.

41

Clearly, acquiring adequate information about the energetics of protein-ion interac-

42

tions is required for an understanding of the biology of Zn metalloproteins. 9 A significant

43

component of the free energy of Zn complexation involves polarization of the ligands

44

and charge transfer. 10,11 Proper estimates of the interaction energies between the ion

45

and the relevant ligands can be derived from quantum mechanical calculations (QM), but

46

high-level QM calculations for a whole protein (or even a protein domain) are still com-

47

putationally prohibited. Approximating some of the system by employing e.g., a quantum

48

mechanics/molecular mechanics (QM/MM) approach (for example, see 12,13 ) can signif-

49

icantly reduce the computational requirements. However, QM/MM calculations present

50

their own challenges, such as the treatment of the boundary between the QM and MM

51

parts, and the computational speed for the QM part, particularly in dynamic systems that

52

contain multivalent metal ions. Moreover, in many such cases the QM region needs to be

53

large and the time-span to be covered need to be relatively long (hundreds of nanosec-

54

onds or longer).

55

A force field (i.e., MM-) based approach is routinely used in molecular dynamics (MD)



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

56

simulations of biological systems. 14,15 Available atomistic, non-polarizable force fields are

57

often fast and efficient but not always accurate enough to represent the intricate character

58

of protein-ion interactions. 16 Some force fields make use of a covalent representation for

59

the zinc-ligand interactions. However, although bonded models can be very useful in the

60

study of certain metalloproteins, 7,17 they are not suited for studies where the metal coor-

61

dination undergoes changes, necessitating the development of non-bonding methods.

62

Several attempts to develop non-bonded force fields for Zn were reported in the lit-

63

erature, 18–21 but modeling different binding conformations and ligand exchange correctly

64

remains a challenge. To exemplify this challenge for catalytic sites in contact with solu-

65

tion, we performed a 5 ns MD simulation of the S100A12 protein with the well established

66

CHARMM27 force field and the zinc ion parameters developed by Stote and Karplus. 18

67

Already after 2 ns of simulation, Zn becomes hydrated, and it adopts an octahedral coor-

68

dination that is not supported by experimental evidence (Figure 1).

69

Polarizable force fields are more complicated and (at least in principle) can be more

70

accurate than classical, non-polarizable ones. Drude-model based force fields that can

71

handle Zn ions have been developed. 22,23 A different approach to representing the polar-

72

izability of atoms is the use of a multipole expansion, as in the AMOEBA (Atomic Multipole

73

Optimized Energetics for Biomolecular Applications) force field. 24 It is important to note,

74

however, that polarizable force field may fail to represent systems in which a significant

75

charge transfer occurs. The so-called Sum of Interactions Between fragments Ab ini-

76

tio computed (SIBFA) 25 approach can deal with charge transfer but is not yet applicable

77

to dynamical metal coordination changes. Another charge transfer polarisable force-field

78

has been developed for simulations of Zn2+ and other ions in water, 26 but does not include

79

parameters for biomolecules, acetate or imidazole.

80

One reason why non-bonded models in force fields tend to fail to correctly represent

81

Zn2+ -protein systems is that the Lennard-Jones (LJ) parameters are usually calibrated

82

for ions in water and not for ions bound to protein-residues. Indeed, accounting for both


Page 4 of 39

Page 5 of 39

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


(a)

(b)

Figure 1: Simulation of S100A12 (PDBID:2WCB 27 ) with CHARMM27 and zinc ion parameters developed by Stote and Karplus; 18 (a) zinc (silver color) in the binding site with three His and one Asp and (b) after 2 ns simulation the zinc ion binding site is clearly more hydrated and the Zn complex adopts an octahedral conformation that is not supported by experiments. 83

ion-water and ion-residue interactions in the parametrisation of the force field was shown

84

to yield better accuracy in simulations of proteins with other ions (Ca2+ , K+ and Na+ ). 28,29

85

Likewise, Sakharov and Lim 11 developed an empirical scheme to handle polarizability and

86

charge transfer for Cys2 His2 Zn binding sites, but this model is not made to be used with

87

other binding conformations.

88

In our previous work, 30 we highlighted the limitations of some widely used QM-methods

89

(in particular DFT functionals) in describing simple, biologically-relevant Zn and Cd-complexes

90

and found that it is desirable to use the second-order Møller-Plesset theory (MP2) method

91

if more accurate alternative such as Coupled Cluster with Single, Double, and Triple ex-

92

citations (CCSD(T)) are not affordable. Here, we examined the variation of interaction

93

energies of structure-optimized Zn-ligand complexes, using MP2, in vacuum and in two

94

solvents described by continuum representations (water and a low dielectric solvent). We

95

compared the interaction energy differences for ten Zn-complexes calculated with QM

96

to those calculated with MM with the aim to elucidate the origin of the differences. The



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

97

amino acid combinations in the complexes represent zinc binding sites from the protein

98

data bank. In addition, we studied complexes where 1-2 ligands were replaced by wa-

99

ter molecules in order to represent partially solvated systems such as expected during

100

catalysis, when the protein binds or releases the Zn2+ ion, or during ligand exchange.

101

The models that were used in this study represent Zn-binding sites of different char-

102

acter with four to six ligands (Figure 2). Complex S1 represents His4 Zn, which has been

103

suggested to be involved in pathogenic Zn2+ catalyzed amyloid aggregation. 31 Model S2

104

represents His3 Glu/Asp as liganding residues. This combination is common among pro-

105

tein structures (e.g., MMP3 32 and S100A12 27 ). We were interested to see the effects of

106

replacing the acetate by a more basic group (methoxide) that represents serine. This is

107

modeled by complex S3. Serine is not a common ligand for Zn in proteins, 6,8 but in some

108

cases it can replace the native ligand without apparent loss in affinity. 33 Complex S4

109

represents His3 Wa (Wa for water). This coordination is found in several enzymes includ-

110

ing carbonic anhydrase 34 and MMP13. 35 Complexes S5–S9 represent intermediates that

111

could be found after ligand exchange, during protein folding and unfolding, or during pep-

112

tide aggregation. Finally, complex S10 is the fully hydrated ion. Cysteine ligands were

113

investigated previously 9,17,36–42 and were therefore not included in this study.

114

Given that many applications and models (e.g. force fields, solvation models) begin

115

with and rely on an accurate description of the gas-phase interaction energies, or at the

116

very least make use of them as a reference, we first performed calculations for nine Zn2+

117

complexes with model amino acid like ligands in vacuum. Thereafter we included also

118

solvent interactions to cover a more complete and realistic description of the interface

119

between solution and protein. 8 In this study we used tetrahydrofuran (THF) as a repre-

120

sentation of the protein environment outside the Zn2+ -ligand complexes, because of the

121

similarity between the dielectric properties of THF and those of proteins. 8,43


Page 6 of 39

Page 7 of 39

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


122

Methods

123

Model systems

124

The Zn-complexes included in this work are shown in Figure 2. The model systems were

125

numbered S1-10, with no particular order. We replaced the amino acid side chains of

126

histidine, serine and glutamic acid (or aspartic acid) by imidazole, methoxide and acetate,

127

respectively. 44,45 Of note, serine was modeled as methoxide and not as methanol because

128

QM calculations have shown that binding of methoxide was favored over water, but that

129

of methanol is not 46 and since metal coordination can reduce the pKa of serine enough

130

to deprotonate the oxygen 47 . Water was also a possible ligand. The carboxylate residue

131

of acetate can interact with the metal ion in a bi-dentate or mono-dentate manner. 48 S10

132

is a reference complex with only water ligands, [Zn(H2 O)6 ]2+

133

Quantum chemistry calculations

134

The GAMESS 49 (General Atomic and Molecular Electronic Structure System) software

135

was used for all QM calculations. Each of the model systems (S1-10) was geometry

136

optimized with the default Quadratic Approximation algorithm with the M06 50 functional

137

and the def2-TZVP 51 basis set. For these optimized structures, interaction energies were

138

then calculated at the MP2 52 level with the aug-cc-pVTZ 53 basis set. Solvent effects were

139

approximated by the iterative Conductor Polarizable Continuum Model (C-PCM). 54 The

140

solvents considered were water (dielectric coefficient r =78.39) and tetrahydrofuran (THF,

141

r =7.58). The atomic radius for Zn was then set to 1.39 Å. 8 For the PCM calculations the

142

aug-cc-pVDZ 53,55 basis set was used.

143

The Counterpoise (CP) correction method 56 was applied for Basis Set Superposi-

144

tion Error (BSSE) correction of interaction energies. The QM interaction energies were

145

further analyzed with the Localized Molecular Orbital Energy Decomposition Analysis



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

Page 8 of 39

S1 Imi4 (+2)

S2 Imi3 Ace (+1)

S3 Imi3 Meo (+1)

S4 Imi3 Wa (+2)

S5 Imi2 Ace (+1)

S6 Imi2 Wa2 (+2)

S7 AceWa2 (+1)

S8 AceWa4 (+1)

S9 AceWa5 (+1)

S10 Wa6 (+2)

Figure 2: Model systems, S1-10. Imi=imidazole (model complex for histidine), Ace=acetate (model complex for aspartic acid or glutamic acid) and Meo=methoxide ion (model complex for serine). Distances ≥ 2.1 Å between the interacting atoms are shown with dashed lines and distances < 2.1 Å are shown as bonds with solid lines. The net charge of the complex is shown in parenthesis.


Page 9 of 39

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


146

(LMOEDA) 57,58 scheme, which divides the interaction energy into contributions from elec-

147

trostatics, exchange, repulsion, polarization (which includes charge transfer) and disper-

148

sion. This Energy Decomposition Analysis (EDA) was carried out in order to shed light

149

on which interactions that play major (or minor) roles in the binding of the metal ion to

150

ligands. The EDAPCM 58 method was used to perform an EDA in solvated systems, using

151

a modified version of GAMESS courtesy of Peifeng Su.

152

Quantum mechanics interaction energies

153

The interaction energy in vacuo between the Zn and all the ligands (L) was calculated for

154

each system (S1-10) according to:

∆E int = Ecomplex − (EZn2+ +

X

EL )

(1)

allL

155

where Ecomplex is the total energy of the complex at its optimized structure. Individual lig-

156

and energies (BSSE-corrected), were calculated at the equilibrium geometry (i.e. lowest

157

energy structure).

158

To explore the influence of ligand-ligand interactions, calculations were also made for

159

the interaction energies of a "ligands-only cluster", i.e. a complex from which the ion has

160

been removed but the ligands were fixed at the positions of the optimized ion-ligands

161

complex, according to:

∆E int,shell = Ecomplex − (EZn2+ + ELC )

(2)

162

Here Ecomplex is the total energy of the complex at its optimized structure, and ELC is the

163

total energy of the "ligands-only cluster" kept at the same geometry as in the optimized Zn-

164

ligand complex. The difference between ∆E int,shell and ∆E int accounts for the influence

165

of ligand-ligand interactions on the formation of the complex.

166

To evaluate the energy cost or gain accompanying the exchange of one or more lig9 ACS Paragon Plus Environment


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

Page 10 of 39

167

ands, the difference between ∆E int,shell for selected complexes were calculated. Such

168

interaction energy differences are given the notation ∆∆E int,shell (SX → SY ) and were

169

calculated as the differences between two selected complexes, SX and SY.

∆∆E int,shell (SX → SY ) = ∆E int,shell (SY ) − ∆E int,shell (SX)

(3)

170

Solvent effects were approximated with the Zn-ligand systems embedded in a polar-

171

izable continuum model (PCM) 54 solvent. Using the PCM model, equation (2) is trans-

172

formed into a free-energy solvent-embedded variant, namely:

CM P CM ∆Gint,shell,P CM = GPcomplex − (GPZnCM ) 2+ + GLC

(4)

173

CM GPcomplex is the total free energy of the complex (using the structure optimized in vacuum)

174

and GPLCCM is the total free energy of the "ligands-only cluster".

175

Optimization of the structures was performed in GAMESS-US with the M06 DFT func-

176

tional. DFT calculations are significantly faster than those performed with MP2, whereas

177

our earlier study has shown that differences with respect to MP2-optimized structures

178

were very small. 30 In the same study, interaction energies calculated with various function-

179

als and with MP2 were compared to CCSD(T) calculations as reference, which revealed

180

good agreement between MP2 and CCSD(T) results (mean absolute error 0.9 kcal/mol),

181

but less convincing agreement between DFT and CCSD(T).

182

Atomic charges

183

We calculated atomic charges by use of three different schemes: Mulliken population

184

analysis, 59 molecular electrostatic potential charge fitting, also known as the Kollman-

185

Singh or Merz-Kollman-Singh (MKS) scheme, 60 and the CHarges from ELectrostatic Po-

186

tentials using a Grid based method (CHELPG). 61 Mulliken’s method is widely used. How-


Page 11 of 39

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


187

ever, its drawbacks include strong dependence on the basis set and its reliance on equal

188

division of overlapping occupied orbitals even between atoms with differences in elec-

189

tronegativities. Electrostatic Potential (ESP) based methods such as MKS and CHELPG

190

do not suffer from these limitations. The main difference between the MKS method and

191

CHELPG is in the procedure used to fit the ESP.

192

Molecular mechanics interaction energy calculations

193

194

Non-polarizable force fields

195

We also performed force-field, i.e. molecular mechanics (MM), calculations of ∆E int,shell

196

for systems S1-S10. Here, the system preparation, energy minimization and interac-

197

tion energy calculations were performed with the GROMACS 62 (GROningen MAchine for

198

Chemical Simulations) software. For the ligands, the SPC/E 63 (extended simple point

199

charge model) water model was used within the OPLS-AA 64,65 (Optimized Potentials for

200

Liquid Simulations-All Atom) and Amber99SB 66 (Assisted Model Building with Energy

201

Refinement) force fields and the TIP3P (three-site Transferable Intermolecular Potential)

202

water model within the CHARMM27 67 force field (see Table 1). For the Zn2+ ion, three

203

different LJ-parameter sets were tested. 11,18,21 We used the combination of force fields

204

and parameters as originally developed, i.e., Li, Merz and co-workers 21 with Amber99SB,

205

Stote and Karplus 18 and Sakharov and Lim 11 with CHARMM27, and we also explored the

206

Zn-parameters from Stote and Karplus with all three force fields. All three LJ-parameter

207

sets were originally optimized for the divalent ion-water systems. Force field parameters

208

for acetate, imidazole and methoxide are presented in Tables S1–S3 in the Supporting

209

information.

210

Coordinates from the QM optimized structures (M06/def2-TZVP) were used as input



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

211

structures for the force field based calculations. The systems were then energy minimized

212

with the conjugate gradient algorithm in vacuum. A plain cutoff at 2.0 nm was used for the

213

Coulomb and van der Waals forces. Table 1: Molecular mechanics set ups used in this study, with non-polarizable force fields

A

LJ-parameter set Reference σZn (nm) Zn (kJ/mol) 18 Stote and Karplus 0.194216 1.0460

B C

Li et al. 21 Sakharov and Lim 11

0.227400 0.156800

0.0148 0.7657

Force field

Water model for ligands OPLS-AA SPC/E Amber99SB SPC/E CHARMM27 TIP3P Amber99SB SPC/E CHARMM27 TIP3P


Page 12 of 39

Page 13 of 39

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


214

215

Polarizable force-field

216

Interaction energies were calculated by applying the CHARMM Drude force field. 22,68 Wa-

217

ter molecules were represented by the polarizable SWM4-NDP model, 69 which is com-

218

patible with the same force field. System preparation was carried out in CHARMM, 70,71

219

based on scripts adapted from CHARMM-GUI. 72,73 Energy minimization and calculations

220

were performed in NAMD. 74,75

221

Results and discussion

222

Interaction energies of Zn complexes - comparison between ab initio

223

and molecular mechanics calculations

224

Complexes S1-S10 include four-, five- and six-coordinated systems (acetate can bind

225

Zn2+ in either a monodentate or bidentate fashion). Optimization of the structures with

226

ab-initio methods resulted in tetrahedral or octahedral coordinations with small distortions

227

(Figure 2). Interestingly, coordination numbers of 5 and 7 were observed for complexes

228

S2 and S9, respectively, with both the Amber99SB and CHARMM27 FFs as both oxygens

229

of the acetate coordinated the Zn2+ ion in the minimum. The conformations of all other

230

systems optimized using the MM methods were very similar to those obtained by DFT.

231

Bond lengths are presented in the Supporting Information Table S4.

232

The interaction energies between the Zn and its ligands (systems S1-S10) as cal-

233

culated with MP2/aug-cc-pVTZ for the DFT-optimized structures are presented in Ta-

234

ble 2. ∆Eint,shell values are presented with CP correction. However, the corrections

235

were 5 kcal/mol or smaller in all systems (Table S5), indicating that the basis set is large

236

enough for these calculations in view of the fact that the interaction energies themselves 13 ACS Paragon Plus Environment


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

237

are large. The least favorable interaction energy (-360 kcal/mol), was calculated for the

238

fully hydrated Zn ion, whereas the most favorable one (-620 kcal/mol) was calculated for

239

complex S3 (Imi3 Meo). Furthermore, a comparison of ∆Eint,shell and ∆E int in Table 2,

240

revealed that the ligand-ligand interactions in the complexes contributed 20 kcal/mol or

241

less to the interaction energies, and were thus fairly insignificant in comparison with the

242

overall magnitude of the interactions. To estimate the influence of the geometric distortion

243

of the ligands in the complex, i.e., deviations from the optimized isolated ligand structures,

244

we recalculated ∆Eint for systems S2, S6, and S10 in the same geometry as in the com-

245

plex rather than at the optimal isolated ligand geometries. The interaction energies were

246

rather similar irrespective of whether optimization of the isolated ligands was performed

247

in all three cases (Table S5).

248

Page 14 of 39

The interaction energies for the ten Zn-ligand systems, calculated with the force fields

249

OPLS-AA, Amber99SB and CHARMM27 in combination with the three different LJ-parameters

250

for Zn, are presented in Table 2 and Figure 3a. The differences between the MM and MP2

251

interaction energies were overall large (65 kcal/mol in average and up to 170 kcal/mol),

252

see Figure 3b. Examination of the three force fields used here reveals that although the

253

actual ∆E int,shell values differed substantially from those calculated with MP2, the three

254

force fields displayed quite similar trends. All overestimated the interaction energies of

255

complexes S8 to S10, and underestimated the rest, where the interaction energies of S1

256

and S4 were the most heavily underestimated. The similarity between the force fields can

257

be attributed to the similarity in partial charges (Table S1-S3 in Supporting Information).

258

Calculations with the Stote and Karplus parameters for Zn in combination with the

259

CHARMM27 force field performed significantly better than the Amber99SB and OPLS-AA

260

force fields for all model systems except S7 (Figure 1). Stote and Karplus developed their

261

parameters to be compatible with the CHARMM series of force fields. As can be seen

262

from Table 2 and Figure 3a, the Zn parameters from Sakharov and Lim, 11 also calculated

263

with the CHARMM27 force field, resulted in better agreement with MP2 than the param-


Page 15 of 39

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


264

eters from Stote and Karplus for complexes S1 and S3-S7, but not for the others. Bond

265

length values are given in the Supporting Information, Table S6. We also note that the

266

optimized Zn parameters from Li et al. (in combination with the Amber99SB force field)

267

led to a better agreement with the MP2 calculation for complexes S1-S7 compared with

268

the LJ-parameters developed by Stote and Karplus, whereas the opposite was true for

269

complexes S8-S10.

270

The calculated interaction energies for complexes S1 to S9, with respect to the hy-

271

drated ion (S10) as a reference, are presented in Table 3. At the MP2 level in vacuum,

272

complexes S1-S9 yielded lower, i.e. more favorable, interaction energies than the hy-

273

drated ion. This was not the case for the interaction energies calculated with the different

274

force fields, where model systems S1, S4 and S6 were less stable than the all-water

275

complex (S10). This may be due to larger contribution from polarization (see section 3.4

276

below).

277

Interaction energies in continuum representations of solvents

278

The QM interaction energies were estimated also in a continuum representation of water

279

and THF (the latter may roughly represent interactions in proteins 8 ). Energy values with a

280

PCM model, formally correspond to free energies of interactions and should therefore not

281

be directly compared with ∆E int,shell . As expected, differences between ∆Gint,shell values

282

in THF and water were seen to be most apparent for systems that include charged ligands

283

(acetate and methoxide), due to the difference in dielectric constant (Figure S2, Support-

284

ing information). For those complexes, the ∆Gint,shell values were about 30 kcal/mol more

285

negative (i.e. more stable) in THF than in water. The differences for the neutral systems,

286

containing only imidazole or water as ligands were much smaller. For the fully hydrated

287

ion, embedding in THF increased the free energy of interaction by as little as 1 kcal/mol,

288

indicating that in PCM is the ion in this sense regarded as fully solvated by a single shell

289

of water molecules. Interestingly, whereas ∆E int,shell values were the least favorable for 15 ACS Paragon Plus Environment


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

290

the hydrated ion (S10, Table 3), this was not the case for ∆Gint,shell . In PCM water, the

291

interaction free energies calculated for complexes S5-S7 were less favorable than S10,

292

which suggests that such conformations would not be stable in solution. One may argue

293

that the low dielectric constant of THF does not well represent the systems that are par-

294

tially hydrated (S6 and S7), and that a higher dielectric constant ( ≥20 or more) would

295

be better. 43 However, we find that these systems yield interaction free energies that are

296

even less favorable in PCM water relative to system S10.

297

Ligand exchange energies

298

To get a measure of the interaction energies for the individual ligands in the complex, the

299

difference between two selected systems, ∆∆E int (SX → SY ) or ∆∆E int,shell (SX → SY )

300

was calculated from the MP2 values. Such ligand exchange interaction energies for (i)

301

exchange from water to imidazole, (ii) from water to acetate, and (iii) from imidazole to

302

acetate are presented in Table 4. An exchange of a water molecule by imidazole, i.e.

303

S6 → S4, is accompanied by an energy gain of about 30–40 kcal/mol, and an exchange

304

of two water molecules by two imidazoles (S7 → S5) by almost twice as much. Ligand

305

exchange of two water molecules by one acetate (S10 → S8) is very favorable in vacuum

306

(-200 kcal/mol), while exchange from two imidazoles to one acetate (S6 → S7) yields

307

-140 kcal/mol.

308

The difference in the interaction energies with regard to ligand exchange in PCM water

309

and THF, ∆∆Gint,P CM (SX → SY ), were similar in magnitude for water to imidazole.

310

However, the exchange from water to acetate yielded only 7–30 kcal/mol in PCM, which

311

was much less than the corresponding value in vacuum. It is apparent that the effects of

312

ligand exchange are less pronounced when a water molecule is replaced by a charged

313

acetate in solvent, especially in water, where the favorable interaction between Zn and the

314

charged ligand is offset by the loss of solvation free energy for the complex who carries out

315

a smaller charge (+1 versus +2 when a water molecule is replaced by an acetate). The 16 ACS Paragon Plus Environment

Page 16 of 39

Page 17 of 39

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

∆∆Gint,shellP CM (SX → SY ) values for a replacement of two imidazoles to one acetate

317

were +40–+70 kcal/mol in PCM, suggesting that such an exchange is not likely to occur.

318

To assess the performance of the three force fields the same ligand exchange ener-

319

gies were calculated by MM. OPLS-AA and Amber99SB yielded no difference in inter-

320

action energies upon exchange from water to imidazole. In contrast, calculations with

321

the CHARMM27 force field resulted in interaction energies which were more favorable

322

by about 20 kcal/mol per molecule. These values are in better agreement with the MP2

323

results (30–40 kcal/mol). It is apparent that the Zn+2 -imidazole (His) interaction is not

324

favorable enough in MM. The distances between the Zn atom and its nitrogen ligands on

325

imidazole are accurate with the OPLS-AA MM, whilst ∼0.10 Å too long with Amber99SB

326

and CHARMM27 (Table S4). Ligand exchange from two water to one acetate is very

327

favorable in vacuum (-200 kcal/mol) also when calculated by MM.

328

For a ligand exchange from two imidazole to one acetate the two force fields (OPLS-

329

AA and Amber99SB) predicted too favorable ∆∆E int,shell (SX → SY ) compared to MP2,

330

whereas the interaction energies calculated with CHARMM27 were closer to MP2.

331

Energy decomposition and atomic charge analyses

332

Energy decomposition analysis 57,58 is used here in order to estimate the relative con-

333

tributions of interactions that govern the binding of Zn in the different complexes. The

334

decomposed interaction energies, are presented in Figure 4 and Table S7.

335

The contribution of polarization to the interaction energies is larger in the systems that

336

involve imidazole (S1-S6) compared to those that do not (S7-S10). The contribution from

337

polarization appears even more significant for S1, S4 and S6 where all ligands are neutral

338

and the total interaction energies are smaller. The electrostatic contribution is roughly as

339

large as all the other stabilizing interactions together for the neutral systems, whereas it

340

is more dominant in all the others. The contribution of dispersion to the total interaction

341

is somewhat more significant for all the imidazole-containing complexes. This appears to 17 ACS Paragon Plus Environment


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

342

play a role in biology since dispersion interactions are relevant for nature’s fine-tuning of

343

biochemical reactions, even if the electrostatic forces are much stronger. 76

344

The interaction energies for S8 and S9 in vacuum are both about -550 kcal/mol, and

345

EDA revealed only minor differences between the systems. The addition of a water ligand

346

is apparently offset by the acetate coordinating in a mono-dentate instead of a bidentate

347

fashion.

348

Turning now to the solvent-embedded (in THF) systems, we find that the electrostatic

349

contributions were larger compared to the calculations in vacuum, whereas the polariza-

350

tion interaction energies became a little bit less significant (Figure 4b). The difference in

351

interaction free energies between the different systems levels off somewhat when they

352

are embedded in a continuum due to the considerably large desolvation energies (Fig-

353

ure S1 in Supporting Information). We note that, the desolvation energy considered by

354

LMOEDA 58 refers only to the electrostatic interaction between solute and solvent in a

355

cavity defined by the complex. Complexes with a negatively charged ligand (S2, S3, S5,

356

S7, S8 and S9) were characterized by desolvation energies of 170–250 kcal/mol. The

357

desolvation energies of complexes with only neutral ligands (S4, S1, S6 and S10) were

358

much lower, between 10 and 80 kcal/mol.

359

The complex with the least favorable interaction energy in PCM is the tetrahedral com-

360

plex S7 with one acetate and two water ligands (Table S8). The high energy cost due

361

to the desolvation for this system resulted in a higher ∆Gint,shell value compared to S8

362

and S9, with their large numbers of water ligands (Figure 4b). It had previously been

363

suggested that a tetrahedral coordination of ligands was preferred for Zn2+ complexes

364

containing one negative ligand. 77 Here we found that hexacoordinated system (S9) had

365

more favorable interaction energy than the tetracoordinated (S7) by 40 kcal/mol in vacuum

366

and by 100 kcal/mol in PCM.

367

Examination of the atomic charge on the Zn ion can yield a qualitative indication of the

368

amount of the electronic charge transfer. In the following, we refer to charges calculated


Page 18 of 39

Page 19 of 39

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


369

by the CHELPG method, which are shown in Table 5 (MKS and Mulliken charges are

370

shown in Table S8). The water molecules apparently polarized the Zn ion, which resulted

371

in a +2.2 au charge for Zn in model system S10, consistent with similar calculations

372

(personal communication with authors of 78 ). In the systems with a single acetate and

373

water ligands (S7-S9) the Zn charge was just over +1 au, whereas in the complexes that

374

contained imidazoles (S1-S6) the charge of the Zn ion ranged between +0.5 and +1 au.

375

The results reveal that charge transfer is evident in all systems. This may indicate that

376

even polarizable force fields are liable to inaccuracies. Inclusions of charge transfer in

377

the calculation of force-field interaction energies 25 may yield better agreement with high-

378

level ab initio results. Several methods have been suggested to overcome this problem

379

by charge fitting. 11,17,79 However, this may not be an adequate solution for systems that

380

contain several His ligands, because (1) the charge transfer was very substantial, (2)

381

there was an evident disagreement even between the two ESP-based methods, and (3)

382

any movement of the ion or ligands is liable to modify the amount of charge transfer.

383

Simulations of the dynamics of such complexes thus remain a challenge.


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

Zn systems Name Ligands S1 Imi4 S2 Imi3 Ace S3 Imi3 Meo S4 Imi3 Wa S5 Imi2 Ace S6 Imi2 Wa2 S7 AceWa2 S8 AceWa4 S9 AceWa5 S10 Wa6

∆Eint -428 -586 -600 -399 -568 -366 -496 -537 -553 -342

∆Eint,shell -438 -606 -619 -409 -571 -376 -513 -555 -554 -357

QM vacuum

QM PCM ∆Gint,shell Water THF -371 -379 -368 -397 -388 -416 -330 -339 -297 -330 -290 -300 -223 -258 -332 -358 -349 -373 -339 -340

MM vacuum (∆Eint,shell ) OPLS-AA Amber99SB CHARMM27 CHARMM-Drude 18 18 21 18 11 A A B A C Polarizable FF -286 -271 -325 -347 -404 -509 -559 -550 -630 -588 -663 -702 -533 -513 -588 -570 -642 N/A -287 -271 -328 -324 -380 -471 -493 -485 -558 -510 -576 -647 -287 -271 -330 -302 -356 -425 -493 -485 -561 -468 -530 -567 -613 -601 -686 -580 -654 -631 -619 -630 -676 -613 -652 -600 -414 -392 -461 -374 -435 -397

Table 2: Interaction energies, ∆Eint and ∆Eint,shell (kcal/mol), calculated with QM (MP2/aug-cc-pTZV) in vacuum and PCM (water and THF), with five combinations of FF and LJ parameters for Zn2+ , see Table 1, and the CHARMM-Drude polarizable force field 22,68 in vacuum for the Zn model systems (S1-S10). N/A - not available.


ACS Paragon Plus Environment

20

Page 20 of 39

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

Zn systems Name Ligands S1 Imi4 S2 Imi3 Ace S3 Imi3 Meo S4 Imi3 Wa S5 Imi2 Ace S6 Imi2 Wa2 S7 AceWa2 S8 AceWa4 S9 AceWa5 S10 Wa6

∆Eint -86 -244 -258 -57 -226 -24 -154 -195 -211 0

∆Eint,shell -80 -249 -262 -52 -214 -19 -156 -198 -197 0

QM vacuum

QM PCM ∆Gint,shell Water THF -32 -39 -29 -57 -49 -76 9 1 42 10 49 40 116 82 7 -18 -10 -33 0 0

MM vacuum (∆Eint,shell ) OPLS-AA Amber99SB CHARMM27 CHARMM-Drude 18 18 21 18 11 A A B A C Polarizable FF 128 121 135 28 31 -112 -145 -158 -169 -214 -228 -305 -119 -121 -127 -195 -207 N/A 127 121 133 50 55 -74 -79 -93 -97 -135 -141 -250 127 121 131 73 79 -28 -79 -93 -100 -93 -95 -170 -199 -209 -225 -206 -218 -234 -205 -238 -215 -239 -218 -203 0 0 0 0 0 0

int,shell int,shell Table 3: Interaction energies for complexes (S1-S9) with respect to hydrated Zn2+ (S10), i.e. ∆ESX − ∆ES10 , 2+ calculated with MP2 in vacuum and in PCM, and with five combinations of FF and LJ parameters for Zn (Table 1) in vacuum. All values are in kcal/mol. N/A - not available.

Page 21 of 39 The Journal of Physical Chemistry

ACS Paragon Plus Environment

21


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

(a)

Page 22 of 39

(b)

Figure 3: Comparison of the MM and QM (MP2) interaction energies (∆Eint,shell ) for complexes S1-S10. The MM energies were calculated with the OPLS-AA, CHARMM27 and Amber99SB for the ligands combined with Zn parameters from Stote and Karplus, 18 Sakaharov and Lim, 11 and Li et al., 21 and with the CHARMM-Drude polarizable force field 22,68 . (a) Correlation plot. (b) Comparison of MM interaction energies relative to QM (MP2) where the reference ("0") is the MP2 energy. All values are in kcal/mol.

384

385

Calculations of Zn-ligand interaction energies with a polarizable force-

386

field

387

Some of the limitations of classical force fields when it comes to interactions between

388

metal ions and proteins may be overcome by the use of polarizable force-field. For this

389

reason, we calculated the interaction energies of complexes S1, S2 and S4–S10 with the

390

CHARMM-Drude polarizable force field. The CHARMM-Drude force field was chosen as

391

it is available in several different MD packages and because it included parameters for

392

all the ligands included in this study except methoxide. The results (Table 2), revealed a

393

marked improvement over all of the non-polarizable force-fields tested here. The MM to

394

QM difference is more consistent compared to the non-polarizable force fields (Figure 3b).

395

Moreover, there is an agreement on the ordering of the interaction energies between the 22 ACS Paragon Plus Environment

Page 23 of 39

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


Table 4: Ligand exchange energies calculated as the difference in interaction energies between two complexes SX and SY, ∆∆E int,shell (SX → SY ), (Eq. 3) denoted diff in the table, calculated with MP2, with OPLS-AA and CHARMM27 force fields with Stote and Karplus parameters for Zn 18 and with a polarizable force field 22,68 in vacuum. The interaction energies between the ion and individual ligands with ligands at their lowest energy structure, ∆∆E int (SX → SY ), (see Eq. 1) for the MP2-method in vacuum were calculated for several complexes and are given in parenthesis. The difference between ∆∆E int (SX → SY ) and ∆∆E int,shell (SX → SY ) is the contribution of ligand-ligand interactions within the complex. MP2 ∆∆Gint,shell (SX → SY ) values are calculated with PCM (water and THF). All energies are in kcal/mol. Wa to Imi

Wa to Imi

S6 Imi2 Wa2 ∆E int,shell -376 (-366) -287 -302 -425

S4 Imi3 Wa ∆E int,shell -409 (-399) -287 -324 -471

diff -33 (-33) 0 -22 -46

MP2 vacuum OPLS-AA CHARMM27 CHARMM-Drude

S7 AceWa2 ∆E int,shell -513 (-496) -493 -468 -567

S5 Imi2 Ace ∆E int,shell -571 (-568) -493 -510 -647

diff -58 (-72) 0 -42 -80



S1 Imi4 ∆E int,shell -438(-428) -286 -347 -509

diff -62 (-62) +1 -45 -84


S10 Wa6 ∆E int,shell -357 (-342) -414 -374 -397


diff -198 (-195) -199 -206 -234



S5 Imi2 Ace ∆E int,shell -571 (-568) -493 -510 -647

diff -195 (-202) -206 -208 -222



diff -137 (-130) -206 -166 -142

MP2 vacuum OPLS-AA CHARMM27 CHARMM-Drude 2 Wa to 2 Imi

2 Wa to Ace

2 Imi to Ace MP2 vacuum OPLS-AA CHARMM27 CHARMM-Drude

S6 Imi2 Wa2 ∆Gint,shell -290 -300

S4 Imi3 Wa ∆Gint,shell -330 -339

diff -40 -39

MP2 PCM water MP2 PCM THF

S7 AceWa2 ∆Gint,shell -223 -258

S5 Imi2 Ace ∆Gint,shell -297 -330

diff -74 -72



S1 Imi4 ∆Gint,shell -371 -379

diff -81 -79


S10 Wa6 ∆Gint,shell -339 -340


diff -7 -18



S5 Imi2 Ace ∆Gint,shell -297 -330

diff -7 -30



diff +67 +42


2 Wa to 2 Imi

2 Wa to Ace

2 Imi to Ace MP2 PCM water MP2 PCM THF



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

S1

S2

211 -219

-84

-72

S5

Electrostatic Exchange

-194 -22

-69

S8

-453

-79

S9

Desolvation

-22

145 -59

-164

-15

a

Dispersion Repulsion

-288

144 -59

-149

Polarization

-186

-478

172

197

S6

-23 -316

S7

-20 -526

203 -81

-208 -92

-218

-18 -481

235

251 -99

-208

-26 -317

S4

S3

181

Page 24 of 39

-168

129

-54

-154

-9

-11 -466

S10

-462


-8 -270

Page 25 of 39

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


252

185 213

204

53

S1

S2

-74

-195 -85

197

S3

-25

-179

-20

-21

-530

-339

-99

-161

-565

234 239

S4

S5

-184

-94

S6

-93

-151

-23

-314

145

251

S8

-59

177

S9

-59

-143

-13 -484

146

194 -119

b

-18

-21

180

-72

-160

-81

-538

-341

S7

203 71

237

65

-9 -484

S10

129 11

-54 -147

-152

-5

-8 -476

-274

Figure 4: Interaction energies for Zn model systems S1-S10 divided into their energy components according to the LMOEDA (in frame a) and EDAPCM (in frame b) schemes with the area of the pie relative to the total interaction energies, a) in vacuum and, b) in PCM THF. Negative values represent a stabilizing interaction, positive values a destabilizing interaction.



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

Table 5: Charges calculated with the CHELPG method at the MP2/aug-cc-pVTZ level for the optimized Zn systems S1-S10. For MKS and Mulliken charges, see Table S8. Zn systems Name Ligands Charge S1 Imi4 +2 S2 Imi3 Ace +1 S3 Imi3 Meo +1 S4 Imi3 Wa +2 S5 Imi2 Ace +1 S6 Imi2 Wa2 +2 S7 AceWa2 +1 S8 AceWa4 +1 S9 AceWa5 +1 S10 Wa6 +2

Zn 0.54 0.77 0.72 0.75 0.98 0.99 1.18 1.20 1.08 2.22

-0.30 -0.33 -0.31 -0.19 -0.32 -0.39 -0.71 -0.74 -0.76 -1.15

CHELPG charges Ligands -0.22 -0.13 -0.06 -0.30 -0.08 -0.76 -0.76 -0.74 -0.21 -0.20 -0.68 -0.24 -0.27 -0.76 -0.35 -0.71 -0.75 -0.27 -0.74 -0.86 -0.69 -0.72 -0.74 -0.75 -0.73 -0.74 -0.75 -0.77 -0.72 -0.78 -0.71 -0.68 -0.69 -0.75 -1.13 -1.11 -1.13 -1.11 -1.15

Color codes: imidazole N, acetate O, methoxide O, water O. 396

MM and QM calculations, which suggests that simulations involving ligand exchange may

397

be possible. Examination of the different complexes reveals that the calculated energies

398

are most similar for complexes S10, S9, S6 and S7 that involve a high degree of hydration

399

(at least two water molecules). Complex S8 is an outlier in this respect - it involves four

400

water molecules and its energy is too favorable by 76 kcal/mol relative to MP2. The rea-

401

son for this is that in the optimized structure the aspartate in complex S8 was bidentate,

402

whereas in S9 it was monodentate and hydrogen-bonded to the water. The bidentate in-

403

teraction was stronger (due to electrostatics) and hence the interaction energy was more

404

favorable for complex S8. Examination of the EDA for complexes S8 and S9 (Figure 4A)

405

reveals that despite the different coordination with respect to the aspartate the contribu-

406

tions of electrostatics, exchange, repulsion, polarization and dispersion are almost the

407

same between those systems. The same is apparently not true in MM, where the two

408

oxygens have the same fixed charge.


Page 26 of 39

Page 27 of 39

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


409

Conclusions

410

In this study of small systems that mimic Zn binding sites in proteins three force fields

411

(OPLS-AA, Amber99SB and CHARMM27) were compared with calculations of the inter-

412

action energies using the MP2 method and a large basis set (which was previously shown

413

to yield good agreement with high level CCSD(T) calculations). In vacuum, all systems

414

with protein ligands (S1-S9) were more stable than the hydrated ion (S10), i.e., they were

415

characterized by more negative interaction energies, using MP2. Complexes with neutral

416

ligands only (imidazole and water) were characterized by less favorable total interaction

417

energies and high amount of polarization according to energy decomposition analysis

418

(LMOEDA), compared to complexes with an acetate ligand. When the systems were em-

419

bedded in PCM the differences in interaction energies between the systems almost dis-

420

appeared due to a large desolvation term for the charged acetate ligand. The Zn-acetate

421

interaction was estimated to be about 200 kcal/mol in vacuum (based on calculations of

422

ligand exchange), whereas it was only 10-30 kcal/mol in water or THF because of the

423

large desolvation term. No such differences were observed for the interaction between

424

Zn and imidazole. Half of the attractive interaction energies in the imidazole containing

425

systems could be attributed to electrostatics and more than one fourth of the interaction

426

to polarization. The interaction free energies for neutral complexes containing imidazole

427

were quite similar regardless of the solvent.

428

The imidazole-containing systems (S1-S6) were characterized by less favorable in-

429

teraction energies in non polarizable force field based (MM) calculations than those cal-

430

culated with MP2, while the three systems without imidazoles result in more favorable

431

interaction energies for MM than MP2. The MM energies are off by at least 20 kcal/mol

432

and up to 150 kcal/mol. Interestingly, the free energy of binding the ion in complexes S1,

433

S4, and S6, that included only neutral ligands, were lower than that of the hydrated ion.

434

Classical, non-polarizable force fields represent the interactions between Zn2+ and its lig-



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

435

ands as a combination of electrostatics and LJ interactions (when using the more flexible

436

non-bonded model). Systems S1, S4, S6 and S10 are those where the electrostatic con-

437

tribution to the binding interaction was lower (Figure 4a). In complexes S1, S4 and S6,

438

the contribution from polarization, exchange and repulsion was rather large. Polarization

439

is only implicitly accounted for by classical force fields (e.g., by tailoring partial charges).

440

Exchange is neglected by force fields, whereas repulsion is represented to some extent

441

by a repulsive term in the LJ interactions, which applies at short interatomic distances.

442

Electrostatics is by far the most dominant contribution in non-polarizable MM force fields,

443

and when it is lower the interactions are less favorable than those calculated by Ab initio

444

methods. Force-field parameters for ions were historically developed for use in aqueous

445

solutions, which is likely to be the cause for the over-stability of highly hydrated complexes

446

(S8, S9 and S10) in the non-polarizable force fields. Unfortunately, a simple shift of the

447

potential so that ion-water interactions are not overly stable is not likely to be the key for

448

simulating Zn2+ dynamically together with peptides while accounting for ligand exchange,

449

and the contribution from exchange and repulsion makes it difficult for polarization alone

450

to fix this issue. It was shown here that the order of the interaction energies is similar

451

when calculated with a polarisable force field and MP2, but that the differences in their

452

magnitude that depend on the system. Ion parameters in the CHARMM-Drude force field

453

were developed by simulating the ions in SWM4-NDP water. 68 One may expect that the

454

generation of a polarizable force-field, e.g. of the Drude type, where special attention

455

during the parametrization procedure is paid to the residues in focus in this study, as well

456

as to methanethiol (for cysteine), would lead to a good agreement between the MM and

457

QM inteaction energies, and thus yield fairly realistic simulations.

458

Supporting Information Available

459

The following files are available free of charge.

Supplementary Information: Comple-


Page 28 of 39

Page 29 of 39

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


460

mentary data including force field input, bond lengths for optimized systems, interaction

461

energies with different basis sets, LMOEDA values and atomic charges with the MKS and

462

Mulliken methods.

463

Acknowledgement

464

The calculations were performed on resources provided by the Swedish National Infras-

465

tructure for Computing (SNIC) at Lunarc, centre for scientific and technical computing for

466

research at Lund University and at PDC center for high performance computing at KTH

467

Royal Institute of Technology (project numbers SNIC 2016/1-55, SNIC 2016/1-222, and

468

SNIC 2016/1-518). This work was supported by the Linnæus Centre of Excellence “Bio-

469

materials Chemistry” (R.F), and by the eSSENCE national strategic research program in

470

e-science (K.H).

471

References

472

473

(1) Maret, W.; Li, Y. Coordination dynamics of zinc in proteins. Chem. Rev. 2009, 109, 4682–4707.

474

(2) Burgess, J. Metal ions in solution; Ellis Horwood: Chichester, 1978.

475

(3) Ohtaki, H.; Radnai, T. Structure and dynamics of hydrated ions. Chemical Reviews

476

477

478

479

480

1993, 93, 1157–1204. ˘ (4) Roe, R. R.; Pang, Y.-P. ZincâAšs Exclusive Tetrahedral Coordination Governed by Its Electronic Structure. Molec. Model. Ann. 1999, 5, 134–140. (5) Harding, M. M. Geometry of metal–ligand interactions in proteins. Acta Cryst. D 2001, 57, 401–411.



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

481

482

(6) Sousa, S. F.; Lopes, A. B.; Fernandes, P. A.; Ramos, M. J. The zinc proteome: A tale of stability and functionality. Dalton Trans. 2009, 38, 7946–7956.

483

(7) Peters, M. B.; Yang, Y.; Wang, B.; Füsti-Molnár, L.; Weaver, M. N.; Merz, K. M. Struc-

484

tural Survey of Zinc Containing Proteins and the Development of the Zinc AMBER

485

Force Field (ZAFF). J. Chem. Theory Comput. 2010, 6, 2935–2947.

486

487

(8) Friedman, R. Structural and computational insights into the versatility of cadmium binding to proteins. Dalton Trans. 2014, 43, 2878–2887.

488

(9) Zhang, J.; Yang, W.; Piquemal, J.-P.; Ren, P. Modeling Structural Coordination and

489

Ligand Binding in Zinc Proteins with a Polarizable Potential. J. Chem. Theory Com-

490

put. 2012, 8, 1314–1324.

491

492

493

494

(10) Hoops, S. C.; Anderson, K. W.; Merz, K. M. Force field design for metalloproteins. J. Am. Chem. Soc. 1991, 113, 8262–8270. (11) Sakharov, D. V.; Lim, C. Zn Protein Simulations Including Charge Transfer and Local Polarization Effects. J. Am. Chem. Soc. 2005, 127, 4921–4929.

495

(12) Corminboeuf, C.; Hu, P.; Tuckerman, M. E.; Zhang, Y. Unexpected Deacetylation

496

Mechanism Suggested by a Density Functional Theory QM/MM Study of Histone-

497

Deacetylase-Like Protein. J Am. Chem. Soc. 2006, 128, 4530–4531.

498

(13) Wu, R.; Hu, P.; Wang, S.; Cao, Z.; Zhang, Y. Flexibility of Catalytic Zinc Coordina-

499

tion in Thermolysin and HDAC8: A Born-Oppenheimer ab Initio QM/MM Molecular

500

Dynamics Study. J. Chem. Theory Comput. 2010, 6, 337–343.

501

502

503

(14) Friedman, R. Aggregation of amyloids in a cellular context: modelling and experiment. Biochem. J. 2011, 438, 415–426. (15) Friedman, R.; Boye, K.; Flatmark, K. Molecular modelling and simulations in cancer


Page 30 of 39

Page 31 of 39

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


504

research. Biochimica et Biophysica Acta (BBA)-Reviews on Cancer 2013, 1836, 1–

505

14.

506

(16) Becconi, O.; Ahlstrand, S.; Salis, A.; Friedman, R. Protein-ion interactions: simula-

507

tions of bovine serum albumin in physiological solutions of NaCl, KCl and LiCl. Isr J

508

Chem 2017, In Press, 0.

509

510

511

512

513

514

(17) Lin, F.; Wang, R. Systematic Derivation of AMBER Force Field Parameters Applicable to Zinc-Containing Systems. J. Chem. Theory Comput. 2010, 6, 1852–1870. (18) Stote, R. H.; Karplus, M. Zinc binding in proteins and solution: a simple but accurate nonbonded representation. Proteins 1995, 23, 12–31. (19) Babu, C. S.; Lim, C. Empirical Force Fields for Biologically Active Divalent Metal Cations in Water. J. Phys. Chem. A 2006, 110, 691–699.

515

(20) Wu, R.; Lu, Z.; Cao, Z.; Zhang, Y. A Transferable Non-bonded Pairwise Force Field

516

to Model Zinc Interactions in Metalloproteins. J. Chem. Theory Comput. 2011, 7,

517

433–443.

518

(21) Li, P.; Roberts, B. P.; Chakravorty, D. K.; Merz, K. M. Rational Design of Particle

519

Mesh Ewald Compatible Lennard-Jones Parameters for +2 Metal Cations in Explicit

520

Solvent. J. Chem. Theory Comput. 2013, 9, 2733–2748.

521

(22) Lopes, P. E.; Huang, J.; Shim, J.; Luo, Y.; Li, H.; Roux, B.; Mackerell, A. D. Force

522

Field for Peptides and Proteins based on the Classical Drude Oscillator. J Chem

523

Theory Comput 2013, 5430–5449.

524

525

(23) Li, P.; Kenneth M. Merz, J. Taking into Account the Ion-Induced Dipole Interaction in the Nonbonded Model of Ions. J. Chem. Theory Comput. 2014, 10, 289–297.

526

(24) Ren, P.; Ponder, J. W. Consistent treatment of inter- and intramolecular polarization

527

in molecular mechanics calculations. J. Comput. Chem. 2002, 23, 1497–1506. 31 ACS Paragon Plus Environment


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

528

(25) Gresh, N.; El Hage, K.; Perahia, D.; Piquemal, J.-P.; Berthomieu, C.; Berthomieu, D.

529

Polarizable molecular mechanics studies of Cu(I)/Zn(II) superoxide dismutase:

530

Bimetallic binding site and structured waters. J. Comput. Chem. 2014, 35, 2096–

531

2106.

532

533

(26) Soniat, M.; Hartman, L.; Rick, S. W. Charge Transfer Models of Zinc and Magnesium in Water. J Chem Theory Comput 2015, 11, 1658–1667.

534

(27) Moroz, O. V.; Blagova, E. V.; Wilkinson, A. J.; Wilson, K. S.; Bronstein, I. B. The

535

Crystal Structures of Human S100A12 in Apo Form and in Complex with Zinc: New

536

Insights into S100A12 Oligomerisation. J. Mol. Biol. 2009, 391, 536 – 551.

537

538

539

540

(28) Project, E.; Nachliel, E.; Gutman, M. Parameterization of Ca+2 -protein interactions for molecular dynamics simulations. J. Comput. Chem. 2008, 29, 1163–1169. (29) Hess, B.; van der Vegt, N. F. Cation specific binding with protein surface charges. Proc. Natl. Acad. Sci. U.S.A. 2009, 106, 13296–13300.

541

(30) Ahlstrand, E.; Spångberg, D.; Hermansson, K.; Friedman, R. Interaction energies

542

between metal ions (Zn2+ and Cd2+ ) and biologically relevant ligands. Int. J. Quant.

543

Chem. 2013, 113, 2554 – 2562.

544

(31) Minicozzi, V.; Stellato, F.; Comai, M.; Serra, M. D.; Potrich, C.; Meyer-Klaucke, W.;

545

Morante, S. Identifying the Minimal Copper- and Zinc-binding Site Sequence in

546

Amyloid-β Peptides. J. Biol. Chem. 2008, 283, 10784–10792.

547

548

549

550

(32) Pelmenschikov, V.; Siegbahn, P. E. M. Catalytic Mechanism of Matrix Metalloproteinases: Two-Layered ONIOM Study. Inorg. Chem. 2002, 41, 5659–5666. (33) Stradal, T. B.; Troxler, H.; Heizmann, C. W.; Gimona, M. Mapping the zinc ligands of S100A2 by site-directed mutagenesis. J. Biol. Chem. 2000, 275, 13219–13227.


Page 32 of 39

Page 33 of 39

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


551

(34) Eriksson, A. E.; Jones, T. A.; Liljas, A. Refined structure of human carbonic an-

552

hydrase II at 2.0 Å resolution. Proteins: Struct., Funct., and Bioinform. 1988, 4,

553

274–282.

554

(35) Engel, C. K.; Pirard, B.; Schimanski, S.; Kirsch, R.; Habermann, J.; Klingler, O.;

555

Schlotte, V.; Weithmann, K. U.; Wendt, K. U. Structural basis for the highly selective

556

inhibition of MMP-13. Chem. Biol. 2005, 12, 181–189.

557

(36) Ryde, U. The coordination chemistry of the structural zinc ion in alcohol dehydroge-

558

nase studied by ab initio quantum chemical calculations. Eur. Biophys. J. 1996, 24,

559

213–221.

560

(37) Tiraboschi, G.; Gresh, N.; Giessner-Prettre, C.; Pedersen, L. G.; Deerfield, D. W.

561

Parallel ab initio and molecular mechanics investigation of polycoordinated Zn(II)

562

complexes with model hard and soft ligands: Variations of binding energy and of

563

its components with number and charges of ligands. J. Comput. Chem. 2000, 21,

564

1011–1039.

565

566

(38) Dudev, T.; Lim, C. Factors Governing the Protonation State of Cysteines in Pro˘ ’ An Ab Initio/CDM Study. J. Am. Chem. Soc. 2002, 124, 6759–6766. teins:âAL

567

(39) De Courcy, B.; Gresh, N.; Piquemal, J.-P. Importance of lone pair interac-

568

tions/redistribution in hard and soft ligands within the active site of alcohol dehydro-

569

genase Zn-metalloenzyme: Insights from electron localization function. Interdiscip.

570

Sci.: Comput. Life Sci. 2009, 1, 55–60.

571

(40) Sakharov, D. V.; Lim, C. Force fields including charge transfer and local polarization

572

effects: Application to proteins containing multi/heavy metal ions. J. Comput. Chem.

573

2009, 30, 191–202.

574

(41) Zhu, X.; Bora, R. P.; Barman, A.; Singh, R.; Prabhakar, R. Dimerization of the Full-



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

575

Length Alzheimer Amyloid β-Peptide (Aβ42) in Explicit Aqueous Solution: A Molec-

576

ular Dynamics Study. J. Phys. Chem. B 2012, 116, 4405–4416.

577

578

579

580

(42) Zhu, T.; Xiao, X.; Ji, C.; Zhang, J. Z. H. A New Quantum Calibrated Force Field for Zinc-Protein Complex. J. Chem. Theory Comput. 2013, 9, 1788–1798. (43) Friedman, R.; Nachliel, E.; Gutman, M. Molecular Dynamics of a Protein Surface: Ion-Residues Interactions. Biophys. J. 2005, 89, 768–781.

581

(44) Cini, R.; Musaev, D. G.; Marzilli, L. G.; Morokuma, K. Molecular orbital study of

582

complexes of zinc(II) with imidazole and water molecules. J. Mol. Struc.: Theochem

583

1997, 392, 55 – 64.

584

(45) Friedman, R.; Nachliel, E.; Gutman, M. Application of Classical Molecular Dynamics

585

for Evaluation of Proton Transfer Mechanism on a Protein. Bioch. Bioph. Acta. 2005,

586

1710, 67–77.

587

(46) Trzaskowski, B.; Adamowicz, L.; Deymier, P. A. A theoretical study of zinc(II) interac-

588

tions with amino acid models and peptide fragments. J. Biol. Inorg. Chem. 2008, 13,

589

133–137.

590

(47) Berthon, G. Critical evaluation of the stability constants of metal complexes of amino

591

acids with polar side chains (Technical Report). Pure & Appl. Chem. 1995, 67, 1117–

592

1240.

593

594

(48) Ryde, U. Carboxylate Binding Modes in Zinc Proteins: A Theoretical Study. Biophys. J. 1999, 77, 2777 – 2787.

595

(49) Schmidt, M. W.; Baldridge, K. K.; Boatz, J. A.; Elbert, S.; Gordon, M.; Jensen, J. H.;

596

Koseki, S.; Matsunaga, N.; Nguyen, K. A.; Su, S. et al. General Atomic and Molecular

597

Electronic Structure System. J. Comput. Chem. 1993, 14, 1347–63.


Page 34 of 39

Page 35 of 39

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

598

599


(50) Zhao, Y.; Truhlar, D. G. Density functionals with broad applicability in chemistry. Acc. Chem. Res. 2008, 41, 157–167.

600

(51) Weigend, F.; Ahlrichs, R. Balanced basis sets of split valence, triple zeta valence and

601

quadruple zeta valence quality for H to Rn: Design and assessment of accuracy.

602

Phys. Chem. Chem. Phys. 2005, 7, 3297–3305.

603

604

(52) Møller, C.; Plesset, M. Note on an approximation treatment for many-electron systems. Phys. Rev. 1934, 46, 0618–0622.

605

(53) Balabanov, N. B.; Peterson, K. A. Systematically convergent basis sets for transition

606

metals. I. All-electron correlation consistent basis sets for the 3d elements Sc–Zn. J.

607

Chem. Phys. 2005, 123, 064107.

608

(54) Li, H.; Jensen, J. H. Improving the efficiency and convergence of geometry opti-

609

mization with the polarizable continuum model: New energy gradients and molecular

610

surface tessellation. J. Comput. Chem. 2004, 25, 1449–1462.

611

612

(55) Feller, D. The role of databases in support of computational chemistry calculations. J. Comput. Chem. 1996, 17, 1571–1586.

613

(56) Boys, S.; Bernardi, F. The calculation of small molecular interactions by the differ-

614

ences of separate total energies. Some procedures with reduced errors. Mol. Phys.

615

1970, 19, 553–566.

616

617

618

619

620

621

(57) Peifeng, S.; Hui, L. Energy decomposition analysis of covalent bonds and intermolecular interactions. J. Chem. Phys. 2009, 131, 014102. (58) Su, P.; Liu, H.; Wu, W. Free energy decomposition analysis of bonding and nonbonding interactions in solution. J. Chem. Phys. 2012, 137, 034111. (59) Mulliken, R. S. Electronic Population Analysis on LCAO-MO Molecular Wave Functions. I. J. Chem. Phys. 1955, 23, 1833–1840. 35 ACS Paragon Plus Environment


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

622

623

(60) Singh, U. C.; Kollman, P. A. An approach to computing electrostatic charges for molecules. J. Comput. Chem. 1984, 5, 129–145.

624

(61) Breneman, C. M.; Wiberg, K. B. Determining atom-centered monopoles from molec-

625

ular electrostatic potentials. The need for high sampling density in formamide con-

626

formational analysis. J. Comput. Chem. 1990, 11, 361–373.

627

628

629

630

(62) Van Der Spoel, D.; Lindahl, E.; Hess, B.; Groenhof, G.; Mark, A. E.; Berendsen, H. J. C. GROMACS: Fast, flexible, and free. J. Comput. Chem. 2005, 26, 1701–1718. (63) Berendsen, H. J. C.; Grigera, J. R.; Straatsma, T. P. The missing term in effective pair potentials. J. Phys. Chem. 1987, 91, 6269–6271.

631

(64) Jorgensen, W. L.; Maxwell, D. S.; Tirado-Rives, J. Development and Testing of the

632

OPLS All-Atom Force Field on Conformational Energetics and Properties of Organic

633

Liquids. J. Am. Chem. Soc. 1996, 118, 11225–11236.

634

(65) Kaminski, G. A.; Friesner, R. A.; Tirado-Rives, J.; Jorgensen, W. L. Evaluation and

635

Reparametrization of the OPLS-AA Force Field for Proteins via Comparison with

636

Accurate Quantum Chemical Calculations on Peptides. J. Phys. Chem. B 2001, 105,

637

6474–6487.

638

(66) Hornak, V.; Abel, R.; Okur, A.; Strockbine, B.; Roitberg, A.; Simmerling, C. Compar-

639

ison of multiple Amber force fields and development of improved protein backbone

640

parameters. Proteins: Structure, Function, and Bioinformatics 2006, 65, 712–725.

641

642

(67) MacKerell, A. D.; Banavali, N.; Foloppe, N. Development and current status of the CHARMM force field for nucleic acids. Biopolymers 2000, 56, 257–265.

643

(68) Yu, H.; Whitfield, T. W.; Harder, E.; Lamoureux, G.; Vorobyov, I.; Anisimov, V. M.;

644

Mackerell, A. D.; Roux, B. Simulating monovalent and divalent ions in aqueous solu-

645

tion using a Drude polarizable force field. J Chem Theory Comput 2010, 6, 774–786. 36 ACS Paragon Plus Environment

Page 36 of 39

Page 37 of 39

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


646

(69) Lamoureux, G.; Harder, E.; Vorobyov, I. V.; Roux, B.; MacKerell Jr, A. D. A polariz-

647

able model of water for molecular dynamics simulations of biomolecules. Chemical

648

Physics Letters 2006, 418, 245 – 249.

649

(70) Brooks, B. R.; Bruccoleri, R. E.; Olafson, B. D.; States, D. J.; Swaminathan, S.;

650

Karplus, M. CHARMM: A Program for Macromolecular Energy, Minimization, and

651

Dynamics Calculations. J. Comp. Chem. 1983, 4, 187–217.

652

(71) Brooks, B. R.; Brooks, C. L.; MacKerell, A. D.; Nilsson, L.; Roux, B.; Won, Y.; Ar-

653

chontis, G.; Bartels, C.; Boresch, S.; Caflisch, A. et al. CHARMM: The Biomolecular

654

Simulation Program. J. Comp. Chem. 2009, 30, 1545–614.

655

656

(72) Jo, S.; Kim, T.; Iyer, V. G.; Im, W. CHARMM-GUI: a web-based graphical user interface for CHARMM. J Comput Chem 2008, 29, 1859–1865.

657

(73) Lee, J.; Cheng, X.; Swails, J. M.; Yeom, M. S.; Eastman, P. K.; Lemkul, J. A.;

658

Wei, S.; Buckner, J.; Jeong, J. C.; Qi, Y. et al. CHARMM-GUI Input Generator for

659

NAMD, GROMACS, AMBER, OpenMM, and CHARMM/OpenMM Simulations Using

660

the CHARMM36 Additive Force Field. J Chem Theory Comput 2016, 12, 405–413.

661

(74) Phillips, J. C.; Braun, R.; Wang, W.; Gumbart, J.; Tajkhorshid, E.; Villa, E.; Chipot, C.;

662

Skeel, R. D.; Kalé, L.; Schulten, K. Scalable molecular dynamics with NAMD. J Com-

663

put Chem 2005, 26, 1781–1802.

664

(75) Jiang, W.; Hardy, D. J.; Phillips, J. C.; Mackerell, A. D.; Schulten, K.; Roux, B. High-

665

performance scalable molecular dynamics simulations of a polarizable force field

666

based on classical Drude oscillators in NAMD. J Phys Chem Lett 2011, 2, 87–92.

667

(76) Tsfadia, Y.; Friedman, R.; Kadmon, J.; Selzer, A.; Nachliel, E.; Gutman, M. Molecular

668

dynamics simulations of palmitate entry into the hydrophobic pocket of the fatty acid

669

binding protein. FEBS Lett. 2007, 581, 1243–1247.



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

670

(77) Dudev, T.; Lim, C. Tetrahedral vs Octahedral Zinc Complexes with Ligands of Biolog-

671

ical Interest: A DFT/CDM Study. J. Am. Chem. Soc. 2000, 122, 11146–11153.

672

(78) Asthagiri, D.; Pratt, L. R.; Paulaitis, M. E.; Rempe, S. B. Hydration Structure and

673

Free Energy of Biomolecularly Specific Aqueous Dications, Including Zn2+ and First

674

Transition Row Metals. J. Am. Chem. Soc. 2004, 126, 1285–1289.

675

(79) Duarte, F.; Bauer, P.; Barrozo, A.; Amrein, B. A.; Purg, M.; Åqvist, J.; Kamerlin, S.

676

C. L. Force Field Independent Metal Parameters Using a Nonbonded Dummy Model.

677

J. Phys. Chem. B 2014, 118, 4351–4362.


Page 38 of 39

Page 39 of 39

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

678


Graphical TOC Entry

679


Interaction Energies in Complexes of Zn and Amino Acids: A

Recommend Documents