Configurational Disorder of Water Hydrogen-Bond Network at the

Jun 21, 2017 - For a more comprehensive list of citations to this article, users are ... Hierarchical lattice models of hydrogen-bond networks in wate...
0 downloads 0 Views 2MB Size
Article The Journal of Physical Configurational Chemistry B is published American Disorder by ofthe Water Chemical Society. 1155 Sixteenth Street N.W., Hydrogen Bond Washington, DC 20036

Subscriber access provided by EAST Published by American TENNESSEE STATE UNIV Chemical Society. Copyright © American Chemical Society.

Network at the Protein Dynamical The Transition Journal of Physical Chemistry B is published

Obaidur Rahaman, Maria by the American Chemical Society. 1155 Kalimeri, Marina Katava, Sixteenth Street N.W., Washington, DC 20036

Subscriber access provided by EAST Published by American TENNESSEE STATE UNIV Chemical Society. Copyright © American Chemical Society.

Alessandro Paciaroni, and Fabio Sterpone The Journal of Physical

J. Phys. Chem. B,Chemistry Just Accepted B is published by the 10.1021/ American Manuscript • DOI: Chemical Society. 1155 acs.jpcb.7b03888 • Publication Sixteenth Street N.W., Date (Web): Washington, 21 Jun 2017 DC 20036

Subscriber access provided by EAST Published by American TENNESSEE STATE UNIV Chemical Society. Copyright © American Chemical Society.

Downloaded from http:// pubs.acs.org on June 23, 2017

Just

The Journal of Physical Chemistry B is published by the American Chemical Society. 1155 Street N.W., Accepted Sixteenth Washington, DC 20036 Subscriber access provided by EAST Published by American TENNESSEE STATE UNIV Chemical Society. Copyright © American Chemical Society.

“Just Accepted” manuscripts have been pe online prior to technical editing, formatting f The Journal of Physical is published as a f Society providesChemistry “Just BAccepted” the American dissemination of by scientific material as soon Chemical Society. 1155 appear in full in PDF format by Sixteenth Streetaccompanied N.W., Washington, DC 20036 fully peer reviewed, but should not be consid Subscriber access provided by EAST Published by American

TENNESSEE STATE UNIV Chemical Society. Copyright © American Chemical Society.

readers and citable by the Digital Object Ide to authors. Therefore, the “Just Accepted” The a Journal of Physical is technic in the journal. After manuscript Chemistry B is published Accepted” Web site and published as an AS by the American Chemical Society. 1155and/or gra changes to the manuscript text Sixteenth Street N.W., Washington, DC 20036 Subscriber access provided by EAST Published by American TENNESSEE STATE UNIV Chemical Society. Copyright © American Chemical Society.

and ethical guidelines that apply to the jo or consequences arising from the use of in The Journal of Physical Chemistry B is published by the American Chemical Society. 1155 Sixteenth Street N.W., Washington, DC 20036 Subscriber access provided by EAST Published by American TENNESSEE STATE UNIV Chemical Society. Copyright © American Chemical Society.

Page 1 The of 21 Journal of Physical Chemistry

1 2 3 4 5 6 7

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

Page 2 of 21

Configurational Disorder of Water Hydrogen Bond Network at the Protein Dynamical Transition Obaidur Rahaman,∗,† Maria Kalimeri,∗,‡ Marina Katava,¶ Alessandro Paciaroni,§ and Fabio Sterpone∗,¶ Institute of Structural Mechanics, Bauhaus-Universitt Weimar, Marienstr. 15, D-99423 Weimar, Germany, Department of Physics, Tampere University of Technology, Tampere, Finland, Laboratoire de Biochimie Th´eorique, IBPC, CNRS UPR9080, Univ. Paris Diderot, Sorbonne Paris Cit´e, 13 rue Pierre et Marie Curie, 75005, Paris, France, and Dipartimento di Fisica e Geologia, Universite di Perugia, via A. Pascoli, 06123 Perugia, Italy E-mail: [email protected]; [email protected]; [email protected]



To whom correspondence should be addressed Institute of Structural Mechanics, Bauhaus-Universitt Weimar, Marienstr. 15, D-99423 Weimar, Germany ‡ Department of Physics, Tampere University of Technology, Tampere, Finland ¶ Laboratoire de Biochimie Th´eorique, IBPC, CNRS UPR9080, Univ. Paris Diderot, Sorbonne Paris Cit´e, 13 rue Pierre et Marie Curie, 75005, Paris, France § Dipartimento di Fisica e Geologia, Universite di Perugia, via A. Pascoli, 06123 Perugia, Italy †

1 ACS Paragon Plus Environment

Page 3 of 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

The Journal of Physical Chemistry

Abstract We introduce a novel strategy to quantify the disorder of extended water-water hydrogenbond (HB) networks sampled in particle-based computer simulations. The method relies on the conformational clustering of the HB connectivity states. We successfully applied it to unveil the fine relationship among the protein dynamical transition in hydrated powder, which marks the activation of protein flexibility at Td ∼ 240K, and the sudden increase of the configurational disorder of the water HB network enveloping the proteins. Our finding links, in the spirit of the Adam-Gibbs relationship, the diffusivity of protein atoms, as quantified by the hydrogen mean square displacements, and the thermodynamic solvent configurational entropy.

Introduction Hydrogen bond (HB) connectivity accounts for water anomalies as pure liquid 1 . Moreover, the structure and dynamics of the HB network controls complex processes that take place in aqueous solution such as the solvation of ions 2 , hydrophobic particles 3 and biomolecules 4 , as well as chemical reactivity 5 . Therefore, characterising the collective behaviour of water HBs and quantifying its configurational states is highly demanded 6–9 , in particular in light of recent experimental results based on advanced spectroscopic techniques 1,10 . A key role of HB network relaxation in giving rise to the so called protein dynamical transition (PDT) has been shown by molecular dynamics (MD) simulations 11 and neutron scattering (NS) measurements 12 . The PDT indicates the sudden non-harmonic increase of the protein atomic fluctuations detected above the critical temperature Td , found in the range of 220 − 240 K. The activation of the hydrogen motions above Td and its non-harmonic character are understood within the framework of a hierarchic multi-minima energy landscape of protein conformations 13 that controls hydrogen jumps in different basins 14 . Experimental evidence suggests a tight connection between the acquired protein flexibility and the response to temperature of the solvating environment 15–17 , 2 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

Page 4 of 21

namely, for a hydrated protein powder, the activation of water translation and rotation 11,12,18 . While the coupling among PDT and the water single particle motion is now well probed, a direct correlation of PDT with the collective reshuffling of the water HB network, and its free energy landscape that ultimately controls this dynamics 19 , is still missing. In this work we present a computational approach to quantify the configurational states of an extended water HB network sampled in particle based computer simulations. The methodology is challenged against the PDT 4,12,14,15,20 . To this purpose we have performed MD simulations of a hydrated lysozyme powder system in the temperature window of 200-300 K. The simulated system contains 8 protein molecules and 2241 water molecules (Nw ), corresponding to a hydration level of 0.35h (h = g water/g Lysozyme), close to the experimental conditions described in Ref. 15 . Our methodology allows to prove the direct correlation among the activation of the protein atomic fluctuations and the water connectivity disorder, a.k.a. the solvent configurational entropy.

Methods System Preparation The powder system containing 8 lysozyme proteins was set up following the procedure designed by Tarek and Tobias 21 . The system was built starting from the crystallographic structure of the Hen Egg Lysozyme protein 2LYM. The 8 proteins were solvated using Nw = 2241 water molecules in order to mimic the experimental condition of h=0.39 for proteins in D2 O, and placed initially in a simulation box of dimensions 76x76x36.5 ˚ A. The protein crystalline structures were disordered by preliminary heating, then the system was equilibrated in the NPT ensemble so as to relax the cell dimensions, the protein packing, and hydration. A pictorial representation of the equilibrated system is given in Figure 1 where the simulation box is replicated using two images along the x and y axis to highlight the molecular packing of the 3 ACS Paragon Plus Environment

Page 5 of 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

The Journal of Physical Chemistry

powder system.

Figure 1: Molecular representation of the simulated system. For sake of clarity the simulation box is replicated using two periodic images along the x and y directions.

MD Simulations The MD simulations were carried out using the NAMD package 22 . The system was modelled using the charmm22/CMAP force field for proteins 23 , and the SPC/E 24 for water that was demonstrated to provide good comparison with Neutron Scattering experiments 21 . The simulations were performed in the NPT ensemble with thermostat τT = 1 ps, and barostat period of 100 f s and decay of 50 f s. Electrostatics was handled by the PME algorithm 25 , with grid space resolution of 1.3 ˚ A and short range cut-off 10 ˚ A. At each temperature, after an equilibration phase of about 10 ns, the production run was extended for extra 10 ns. In order to control the convergence of our results, for a few selected temperatures (200, 240, 260 K) we have extended the simulations to 100 ns. 4 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

Page 6 of 21

Analysis The water dynamics was monitored by computing the rotational characteristic time hτ2 i extracted from the second order Legendre polynomial time correlation function for the water dipole µ, P2 = 12 h3 cos2 θ(t) − 1i, where cosθ(t) = µ(t) · µ(0). For convenience, the P2 (t) was fitted using two exponential functions associated to two different time scales, tens (τα ), and hundreds (τβ ) of picoseconds, and a constant parameter accounting for the contribution of immobilized water at the simulation timescale, see SI. The hτ2 i is obtained by averaging the first two characteristic times that are pertinent to the timescale accessible by the NS experiments. The water HB connectivity was analyzed via a clustering procedure performed with in house developed codes. The clustering exploit the partitional algorithm described in Ref. 26 . The details of the methods are described in the Results section. The code is available at https://github.com/mkalime/GromosNetworkClustering.

Results Protein Dynamical Transition In Fig.

2 we report the temperature evolution of the mean square displacement (MSD)

computed for protein hydrogen atoms attached to carbons and averaged over the eight proteins. For each protein, the MSD was calculated directly from the trajectories as hu2 i = P h N1 i |ri (t + τ ) − ri (t)|2 i, where the index i runs over the selected hydrogens and the brackets indicate time average. The time window τ = 150 ps, corresponding to an energy resolution of ∼9 µeV in the NS experiment 15 , was used for computing the displacement. The results from simulations agree with the experimental data obtained by Paciaroni et al. 15 for hydration level h=0.4, with a break of the hydrogen MSD occurring between 230 and 250 K, which is the signature of the PDT 11,12,18,27,28 . This break is clearly visible by inspecting the derivative

5 ACS Paragon Plus Environment

Page 7 of 21

-1

1/ [ps ]

2

1.5

2

MSD (A )

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

1

0.03 0.025 0.02 0.015 0.01 0.005 0 200 220 240 260 280 300 T [K] τ=150 ps exp h=0.3 exp h=0.4

0.5

0 150

200

250

300

T (K)

Figure 2: Mean square displacement (MSD) computed for hydrogens attached to carbon atoms as a function of temperature. NS experimental data for lysozyme powder system at hydration level h = 0.3 and h = 0.4 are shown as stars 15 . The line is a guide to the eye. In the inset plot we report the temperature variation of the average reorientational time hτ2 i, that also shows a crossover at about 240K. of the MSD with respect to the temperature as it is reported in Fig. S1 in SI. The standard deviation associated to MSD is in the order of 3%. Our results are well converged as proved by extending for selected temperatures the trajectories up to 100 ns, and by computing the MSDs considering independent blocks of 10 ns each, see Table S1 in SI. In Fig. S2 we also report the difference between the protein MSDs at 260K and at 230K, i.e. across the dynamical transition, for all the proteins in the system along the respective aminoacidic sequence. It is clear that the local response to the thermal excitations depends on the chosen protein and is a function of the packing in the powder matrix. In agreement with literature 11 , the inset of Fig. 1 shows that this protein dynamic activation takes place at the onset of the reorientation motions of water molecules in the hydration shell, as probed by characteristic reorientation time hτ2 i, see Methods. The activation of water motion around Td ' 240 K suggests an increasing disorder of the HB water layers surrounding the proteins. However, the thermally activated single-particle water reorientation, that can be fruitfully described in the framework of the extended jump model for HB switch events 29,30 , provides only a local kinetic picture. A thorough approach to 6 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

Page 8 of 21

characterise the reorganization of the collective structure of the water-water HB network that extends through the powder matrix would reinforce our understanding.

Disorder of HB connectivity Water is a connected liquid living above the percolation threshold. The response of the HB connectivity to boundaries as well as thermodynamic conditions alters this percolation propensity 4 . For example, the percolation is lost within the hydration layer of an extended rough hydrophobic surface 31 . Similarly, by raising the temperature above the ambient condition, ∼ 310 K, the HB network of water in the hydration layer of small peptides disconnects 32 . In order to follow this percolation propensity in the hydrated powder we have constructed the Nw∗ × Nw∗ connectivity matrix Ai,j (t) for the water molecules surrounding the proteins, Nw∗ (almost 97% of all water in the system). We defined the hydration shell by considering a spherical cut-off of 4.5 ˚ A from the protein heavy atoms. The matrix element Ai,j (t) at the time t takes the value of 1 if molecules i and j are hydrogen bonded, and zero otherwise. The connectivity is established by using the standard geometrical criteria for water HB, i.e. the d < 30o 33 . From the connectivity oxygen-oxygen distance less than 3.5 ˚ A and the angle OOH matrix we calculated the size Smax of the largest cluster formed by all HB connected molecules. In Fig. 3 we report the probability distribution of the percolation fraction p = Smax /Nw∗ . For a percolating structure, the parameter p approaches the upper limit value of 1, on the other hand, it approaches zero when all entities are disconnected. Quite strikingly the percolation fraction responds to temperature in the same way as the single particle mobility. Indeed, for T < 240 K all the probability distributions of p are centered at value p ∼ 0.91. Above the critical temperature 240 K the distribution widens and gradually shifts toward a smaller value of p. This little but very clear relative change of p values distributions (centered from 0.91 to 0.88), indicates a clear activation of the disorder in the water HB connectivity.

7 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

The Journal of Physical Chemistry

probability

Page 9 of 21

0.84

0.86

0.88

0.90

0.92

0.94

0.96

0.98

1.00

percolation fraction

Figure 3: Probability distribution of the percolation propensity of the hydration water at different temperatures.

The configurational landscape of water HB network The disorder of the HB network can be better quantified by exploiting concepts from graph theory 34–36 . A given HB network can be viewed as a dynamic graph G(t), where the water molecules are the vertexes V (G), and the connecting HBs its edges E(G). The graph G is timedependent as the HB network around the proteins is a dynamic entity with a variable number of vertexes and edges. However, for each temperature, it maintains some core topological characteristics due to the condensed phase nature of the liquid water. In order to explore the connectivity states sampled by water in our system, at each temperature we group together the different configurations of the HB network according to structural similarity. For that purpose we employ a clustering analysis that is adapted to the case of dynamic graphs with variable number of vertexes, and which follows in spirit the work of Ref. 35 . Characterisation of the network of HBs in water or formed by ions in an aqueous solution was also recently proposed using graph theory concepts 36,37 without however attempting clustering. As a first step we introduce a convenient order parameter that characterises the structure of the network. This is the eigenvalue spectrum of the Laplacian matrix of the graph G 34,35 . For a given time t, the Laplacian is a Nw∗ × Nw∗ matrix L := L(G) = (Lij ), i, j = 1, ..., Nw∗ , obtained

8 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

Page 10 of 21

from the connectivity matrix Aij , or adjacency matrix in graph theory language, and reads as:

Lij =

    1    

−       0

if ni = nj , √

1 dni dnj

if ni nj ∈ E(G),

(1)

otherwise.

where dni is the number of edges linking ni , for any vertex ni ∈ V (G). In the expression of the Laplacian the HB connecting two molecules (the vertexes ni and nj ) is weighted by the surrounding connectivity, i.e. the number of extra links formed by the two molecules. Therefore, the Laplacian contains in a compact way the information of the global structure of the network, and its eigenvalue spectrum is a suitable order parameter to describe the water connectivity 34 . In fact, the spectrum of the symmetric normalised Laplacian has been shown to capture very well the structural properties of a network 35 and in cases even better than the spectrum of other representations, such as the adjacency matrix which is influenced by the degree distribution of the network 38,39 . Furthermore, its mathematical property to have eigenvalues in the close interval [0-2] allows for a direct pairwise comparison between different networks in time, where even the number of nodes (water molecules) may be variable. Therefore, this allowed us to compare the time evolution of the water connectivity in a finite volume. A probability density estimate for the spectrum is built using a gaussian kernel, and further used to define a metric to compare different network structures. For sake of exemplarity we have plotted in top of Fig. 4 the Laplacian eigenvalue spectra for two water HB connectivities in the hydration shell of a target amino acid of the protein,the Ala42 of one of the proteins in the system. This residue is well exposed to water allowing for a good characterization of its hydration shell. We have individuated two extreme connectivities: in the first one almost all water molecules are connected forming a large cluster of 9 molecules, an isolated HB pair and a disconnected water (see red label on the top-left of the figure), in the

9 ACS Paragon Plus Environment

Page 11 of 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

The Journal of Physical Chemistry

second one, two disconnected large clusters are formed by 4 and 5 water molecules each, while there are additionally 3 isolated molecules and a single HB pair (see yellow label on the topright of the figure). The associated spectra show two typical shapes, a bimodal distribution that characterises a large connected structure (red) and a distribution peaked around the 1 eigenvalue, that is typical of disconnected structure (yellow). The molecular structures of these connectivities are also reported in the corners of the top part of the figures. The interested reader can find a discussion on the relationship among the behaviour of the eigenvalue spectrum and the geometrical features of the network in Ref. 35 , also in SI we provide for a set of ideal example connectivities the associated kernel densities of the Laplacian spectra that can be used as fingerprints to reconstruct the behaviour of more complex cases, see Fig. S3. Additionally, the time evolution of the Laplacian eigenvalue spectrum can be followed so as to obtain information about the characteristic time of HB network structural reorganization, see bottom chart in Fig. 4. In the time evolution we observe the spectrum fluctuating around a structure weakly peaked at the eigenvalue of one with the presence of bimodal shoulders, the spectrum takes at 4 ps the pure bimodal shape characteristic of a largely connected structure, and suddenly after, at 6 ps, shows the transition to a strongly peaked shape characteristic of very disconnected HB structures. This transition clearly encodes HB switch events that lead to the reshuffling of the overall connectivity as well as the penetration of molecules in the hydration shells due to water dynamics or the amino-acid motion. Subsequently we used the Jensen-Shannon divergence 40 to pairwise compare the probability densities obtained from the eigenvalue spectra. For two probability distributions p1 and p2 , associated to two HB network configurations, the Jensen-Shannon divergence is defined as: 1 1 J(p1 , p2 ) = I(p1 , p) + I(p2 , p), 2 2

10 ACS Paragon Plus Environment

(2)

0.4 0.3 0.2 0.0

0.1

density

Page 12 of 21

1(3)-2-4-5

1-2-9

0.0

0.5

1.0

1.5

2.0

eigenvalue

0.7

probability

0.6 0.5 0.4

it

0.3 0.2 0.1

d

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

0.5

The Journal of Physical Chemistry

0.0

Figure 4: Water-Water HB networks enveloping the residue Ala42 of chain A. Top. Two representative states of the HB connectivity are considered. The associated spectral eigenvalues densities of the Laplacian matrix are plotted in the main chart for both connectivities. A molecular view of the hydration shell is provided as inset along with the size of the detected HB clusters, e.g. for the case indicated in red, in the hydration shell we find at a given time one disconnected water, a dimer and a cluster of 9 molecules. Bottom: Spectral diffusion for a time window of 10 ps. where p = 12 (p1 + p2 ) and I(p1(2) , p) =

X

p1(2) (x) log

x∈X

p1(2) (x) . p(x)

(3)

Having a suitable metric in hand, we use a partitional algorithm to cluster the network configurations on the basis of a numerical cut-off on the distance J 26 , and estimate the number of HB network configurations generated by the water molecules at the different temperatures. We quantify this via the total number of clusters found for each temperature, which is shown in panel a of Fig. 5 for a given cut-off Jc = 0.1. This number displays a sudden increase manifested in a peak at 240 K, sensing, once again, the temperature of the protein dynamical transition Td . The growth of the number of clusters signifies that, above the transition temperature,

11 ACS Paragon Plus Environment

Page 13 of 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

The Journal of Physical Chemistry

the HB network begins sampling a larger number of connectivity states. Therefore the protein dynamical transition is directly connected to a sudden increase of the configurational entropy of the water. A graph representation of the configurational landscape for the water HB network is provided in panel c in Fig. 5 for selected temperatures, from which one can appreciate the change in the number of states and their interconnections. On the basis of the transition matrix obtained describing the interconversion from different connectivity states it is possible to perform markovian kinetic clustering 41 so as to cast fast configurational changes and to single out the structural states of the HB network that are separated by high free energy barrier. We provide an example in SI. We have considered also a more strict definition of the HB using d < 20o ), different parameters for the O−O distance (i.e. 3.2 ˚ A) and the HB directionality (OOH and obtained qualitatively the same results. From the network it is possible to estimate the normalised Shannon entropy associated P to the configurational changes of the HB connectivity, Sc = − i xi lnxi , with xi the relative occupancy of the cluster i in the network, see inset of panel a Fig. 5. It is worth to note the deviation from the monotonic behaviour represented by the presence of a sharp peak of Sc at the PDT temperature, between temperatures 220 and 240 K corresponding to ∆Sc /kb ≈ 0.13, and due to the larger number of connectivity states. For a more strict geometrical definition of the water-water HB, see above, the entropy difference is very similar, ∆Sc /kb ≈ 0.10. In the case of supercooled pure water, the growth of configurational entropy was reported to be monotonic without showing any break as a function of temperature 19 . Thus, the sudden increase of entropy at Td , and the associated peak, are clearly related to interfacial nature of water in the protein powder and its coupling with the biomolecular surface. Moreover, the jump in Sc corresponds to the reported peak in the specific heat Cp calculated for hydrated proteins and DNA, and ascribed as well to the coupling between the solvent and the biomolecule 42,43 . Whether this peak relates to a true phase transition of the confined liquid is still debated 42 . A fit according to Sc = K(1/TK −1/T ) 44 , see inset of panel a of Fig. 5, allows to estimate the

12 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

Page 14 of 21

Figure 5: Top panel. Left chart: Number of configurational states of the water-water HB network as a function of temperature. The clustering procedure was performed using a cut-off value Jc = 0.1. Right chart: Adam-Gibbs relationship for the single particle reorietational time and the configurational entropy Sc estimated for the HB network states. Bottom panel: A representation of the ensemble of clusters is given for selected temperatures in the form of network. The size of each vertex is proportional to the configurational state occupancy while the edges represent configurational states interconversions. Kauzmann temperature where the water configurational entropy vanishes, TK = 153 K, very close to what was reported for pure water by Stanley and coworkers 19 . Note that although the absolute value of the number of states change, the observed trend in temperature is independent of the choice of the cut-off value Jc used for clustering if the appropriate range of values is chosen, see Fig S4 in SI. In fact for values of Jc ∈ [0.05, 0.15] the entropy variation is very similar. For a too small value Jc = 0.025 too many clusters get individuated loosing meaningful structural information as effect of statistical noise. A too high threshold generates a few number of clusters by merging together very different connectivity structures. From the physical point of view a fine tune of the parameter Jc can be done by comparing with the configurational entropy estimated by other thermodynamics methods as in Ref. 19 . Beside this technical aspect, the clustering procedure allows also to link, as a function of temperature, the water HB configurational entropy and the single particle kinetics. Very inter-

13 ACS Paragon Plus Environment

Page 15 of 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

The Journal of Physical Chemistry

estingly, above the Td , we found an Adams-Gibbs-like exponential relationship among the water reorientational times and the configurational entropy Sc , hτ2 i−1 ∝ e[A/T Sc ] , similarly to what previously found valid for bulk supercooled water by considering translational diffusion 19,45 , see panel b of Fig. 5, with A a constant that is proportional to the free energy barrier for the network reorganization and is related to the intermolecular potential energy. We consider a final aspect concerning the putative coupling between the protein dynamical transition and a clear transition in water density as it was put forward by analysing NS data 46 , although criticised later by others 47 . In our simulations, the reconfiguration of the water HB connectivity above Td comes along with a decrease of the water density quantified by Voronoi tessellation 48 , but we do not see any signature of the coexistence of separate high and low density populations as the distribution of the single molecule volume is unimodal, see Fig. S6, as previously found in bulk water 49 .

Conclusion To conclude, we have devised a novel procedure to cluster the connectivity states of the water HB network sampled in particle-based computer simulations that allows to access the configuration entropy of the liquid. This allowed us to highlight an intimate correlation among the protein dynamical transition and the increase of the configurational entropy of the water HB connectivity. We show that hydration water acts as an entropic reservoir to sustain protein local conformational changes 50 that play a key role in functionality. Our findings reinterpret the connection between PDT and the activation of water reorientation and translation 12,18 in terms of the Adam-Gibbs microscopic picture 19,45 . Apart from providing a precious piece of information about the relationship between kinetics and thermodynamics at the protein/solvent interface, our results support the view of a connection with the collective structural relaxation of the solvent 16 . The graph theory-based approach successfully applied here could be a key tool

14 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

Page 16 of 21

to study the role of solvent connectivity in complex processes, like ice nucleation, long range proton transfer, the impact of glassy matrices on biomolecules dynamics and stability. It can also be deployed to investigate the clustering of colloidal particles or proteins in suspension.

Supporting Information SI includes data for convergence of MSD, fit parameters for water rotational dynamics, data for HB network clustering, volume analysis.

Acknowledgement The research leading to these results has received funding from the ERC (FP7/2007-2013) Grant Agreement no.258748. Part of this work was performed using HPC resources from GENCI [CINES and TGCC] (Grant x201676818). We acknowledge the financial support for infrastructures from ANR-11-LABX-0011-01 (Dynamo).

References (1) Pettersson, L. G. M.; Henchman, R. H.; Nilsson, A. Water—The Most Anomalous Liquid. Chem. Rev. 2016, 116, 7459–7462. (2) Marcus, Y. Effect of Ions on the Structure of Water: Structure Making and Breaking. Chem. Rev. 2009, 109, 1346–1370. (3) Chandler, D. Interfaces and the Driving Force of Hydrophobic Assembly. Nature 2005, 437, 640–647. (4) Bellissent-Funel, M.-C.; Hassanali, A.; Havenith, M.; Henchman, R.; Pohl, P.; Ster-

15 ACS Paragon Plus Environment

Page 17 of 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

The Journal of Physical Chemistry

pone, F.; van der Spoel, D.; Xu, Y.; Garcia, A. E. Water Determines the Structure and Dynamics of Proteins. Chem. Rev. 2016, 116, 7673–7697. (5) Hynes, J. T. Molecules in Motion: Chemical Reaction and Allied Dynamics in Solution and Elsewhere. Ann. Rev. Phys. Chem. 2015, 66, 1–20. (6) Rao, F.; Garrett-Roe, S.; Hamm, P. Structural Inhomogeneity of Water by Complex Network Analysis. J. Phys. Chem. B 2010, 114, 15598–15604. (7) Heyden, M.; Sun, J.; Funkner, S.; Mathias, G.; Forbert, H.; Havenith, M.; Marx, D. Dissecting the THz Spectrum of Liquid Water from First Principles via Correlations in Time and Space. Proc. Natl. Acad. Sci. USA 2010, 107, 12068–12073. (8) Martin, D. R.; Matyushov, D. V. Dipolar Nanodomains in Protein Hydration Shells. J. Phys. Chem. Lett. 2015, 6, 407–412. (9) Hamm, P. Markov State Model of the Two-State Behaviour of Water. J. Chem. Phys. 2016, 145, 134501. (10) Perakis, F.; Marco, L. D.; Shalit, A.; Tang, F.; Kann, Z.; K¨ uhne, T. D.; Torre, R.; Bonn, M.; Nagata, Y. Vibrational Spectroscopy and Dynamics of Water. Chem. Rev. 2016, 116, 7590–7607. (11) Tarek, M.; Tobias, D. J. Role of Protein-Water Hydrogen Bond Dynamics in the Protein Dynamical Transition. Phys. Rev. Lett. 2002, 88, 138101. (12) Schiro, G.; Fichou, Y.; Gallat, F.-X.; Wood, K.; Gabel, F.; Moulin, M.; H¨artlein, M.; Heyden, M.; Colletier, J.-P.; Orecchini, A. et al. Translational Diffusion of Hydration Water Correlates with Functional Motions in Folded and Intrinsically Disordered Proteins. Nat. Comm. 2015, 6, 6490.

16 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

Page 18 of 21

(13) Frauenfelder, H.; Sligar, S.; Wolynes, P. The Energy Landscapes and Motions of Proteins. Science 1991, 254, 1598–1603. (14) Doster, W.; Cusack, S.; Petry, W. Dynamical Transition of Myoglobin Revealed by Inelastic Neutron Scattering. Nature 1989, 337, 754–756. (15) Paciaroni, A.; Cinelli, S.; Cornicchi, E.; Francesco, A. D.; Onori, G. Fast Fluctuations in Protein Powders: The Role of Hydration. Chem. Phys. Lett. 2005, 410, 400–403. (16) Cornicchi, E.; Onori, G.; Paciaroni, A. Picosecond-Time-Scale Fluctuations of Proteins in Glassy Matrices: The Role of Viscosity. Phys. Rev. Lett. 2005, 95, 158104. (17) Capaccioli, S.; Ngai, K.; Ancherbak, S.; A. Paciaroni, A. Evidence of Coexistence of Change of Caged Dynamics at Tg. J. Phys. Chem. B 2012, 116, 1745–1757. (18) Wood, K.; Frolich, A.; Paciaroni, A.; Moulin, M.; Hartlein, M.; Zaccai, G.; Tobias, D. J.; Weik, M. Coincidence of Dynamical Transitions in a Soluble Protein and Its Hydration Water: Direct Measurements by Neutron Scattering and MD Simulations. J. Am. Chem. Soc. 2008, 130, 4586–4587. (19) Scala, A.; Starr, F. W.; Nave, E. L.; Sciortino, F.; Stanley, H. Configurational Entropy and Diffusivity of Supercooled Water. Science 2000, 406, 166–169. (20) Engler, N.; Ostermann, A.; Niimura, N.; Parak, F. G. Hydrogen Atoms in Proteins: Positions and Dynamics. Proc. Natl. Acad. Sci. USA 2003, 100, 10243–10248. (21) Tarek, M.; Tobias, D. J. The Dynamics of Protein Hydration Water: A Quantitative Comparison of Molecular Dynamics Simulations and Neutron-scattering Experiments. Biophys. J. 2000, 79, 3244–3257. (22) Phillips, J. C.; Braunand, R.; Wang, W.; Gumbart, J.; Tajkhorshid, E.; Villa, E.;

17 ACS Paragon Plus Environment

Page 19 of 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

The Journal of Physical Chemistry

Chipot, C.; ; Skeel, R. D.; Kal´e, L. et al. Scalable Molecular Dynamics with NAMD. J. Comp. Chem. 2005, 26, 1781–1802. (23) MacKerell, A. D.; Feig, M.; Brooks(III), C. L. Extending the Treatment of Backbone Energetics in Protein Force Fields: Limitations of Gas-phase Quantum Mechanics in Reproducing Protein Conformational Distributions in Molecular Dynamics Simulations. J. Comput. Chem. 2004, 25, 1400–1415. (24) Berendsen, H. J. C.; Grigera, J. R.; Straatsma, T. P. The Missing Term in Effective Pair Potentials. J. Phys. Chem. 1987, 91, 6269–6271. (25) Darden, T.; York, D.; Pedersen, L. Particle Mesh Ewald: An Nlog(N) Method for Ewald Sums in Large Systems. J. Chem. Phys. 1993, 98, 10089–10092. (26) Daura, X.; Gademann, K.; Jaun, B.; Seebach, D.; Van Gunsteren, W.; Mark, A. Peptide Folding: When Simulation Meets Experiment. Angew. Chem. Int. Ed. (English) 1999, 38, 236–240. (27) He, Y.; Ku, P. I.; Knab, J. R.; Chen, J.; Markelz, A. G. Protein Dynamical Transition Does Not Require Protein Structure. Phys. Rev. Lett. 2008, 101, 178103. (28) Kim, S. B.; Gupta, D. R.; Debenedetti, P. G. Computational Investigation of Dynamical Transitions in Trp-Cage Miniprotein Powders. Sci. Rep. 2016, 6, 25612. (29) Laage, D.; Hynes, J. A Molecular Jump Mechanism of Water Reorientation. Science 2006, 311, 832–835. (30) Laage, D.; Stirnemann, G.; Sterpone, F.; Rey, R.; Hynes, J. Reorientation and Allied Dynamics in Water and Aqueous Solutions. Ann. Rev. Phys. Chem. 2011, 62, 395–416. (31) Pizzitutti, F.; Marchi, M.; Sterpone, F.; Rossky, P. J. How Protein Surfaces Induce Anomalous Dynamics of Hydration Water. J. Phys. Chem. B 2007, 111, 7584–7590. 18 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

Page 20 of 21

(32) Oleinikova, A.; Brovchenko, I. What Determines the Thermal Stability of the HydrogenBonded Water Network Enveloping Peptides? J. Phys. Chem. Lett. 2011, 2, 765–769. (33) Ferrario, M.; Haughley, M.; McDonald, I. R.; Klein, M. L. MolecularDynamics Simulation of Aqueous Mixtures: Methanol, Acetone, and Ammonia. J. Chem. Phys. 1990, 93, 5156– 5166. (34) Chung, F. R. K. Spectral Graph Theory; Amer. Math. Soc.: Providence, RI, 1997. (35) Banerjee, A. Structural Distance and Evolutionary Relationship of Networks. BioSystems 2012, 107, 186–196. (36) Bak´o, I.; Bencsura, A.; Hermannson, K.; B´alint, S.; Gr´osz, T.; Chihaia, V.; Ol´ah, J. Hydrogen Bond Network Topology in Liquid Water and Methanol: a Graph Theory Approach. Phys. Chem. Chem. Phys. 2013, 15, 15163–15171. (37) Choi, J.-H.; Cho, M. Ion Aggregation in High Salt Solutions. II. Spectral Graph Analysis of Water Hydrogen-Bonding Network and Ion Aggregate Structures. J. Chem. Phys. 2014, 141, 154502. (38) Banerjee, A. The Spectrum of the Graph Laplacian as a Tool for Analyzing Structure and Evolution of Networks. M.Sc. thesis, PhD Thesis. Von der Fakultat fur Mathematik und Informatik der Universitat Leipzig, 2008. (39) Zhan, C.; Chen, G.; Yeug, L. On the Distributions of Laplacian Eigenvalues Versus Node Degrees in Complex Networks. Physica A 2010, 389, 1779–1788. (40) Lin, J. Divergence Measures Based on the Shannon Entropy. IEEE Transactions on Information Theory 1991, 37, 145–151. (41) van Dongen, S. M. Graph Clustering by Flow Simulation. Ph.D. thesis, University of Utrecht, The Netherlands, 2000. 19 ACS Paragon Plus Environment

Page 21 of 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

The Journal of Physical Chemistry

(42) Kumar, P.; Xu, Y. L.; Mazza, M. G.; Buldyrev, S.; Chen, S.-H.; Sastry, S.; Stanley, H. E. Glass Transition in Biomolecules and the Liquid-Liquid Critical Point of Water. Phys. Rev. Lett. 2006, 97, 177802. (43) Mazza, M. G.; Stokelya, K.; Pagnotta, S. E.; Bruni, F.; Stanley, H. E.; Franzese, G. More than one Dynamic Crossover in Protein Hydration Water. Proc. Natl. Acad. Sci. USA 2011, 108, 19873–19878. (44) Angell, C. A. Entropy and Fragility in Supercooling Liquids. J. Res. Natl. Inst. Stand. Technol. 1997, 102, 171–185. (45) Adams, G.; Gibbs, J. H. On the Temperature Dependence of Copperative Relaxation Properties in Glass-Forming Liquids. J. Chem. Phys. 1958, 43, 139–146. (46) Chen, S.-H.; Liu, L.; Fratini, E.; Baglioni, P.; Faraone, A.; Mamontov, E. Observation of Fragile-to-Strong Dynamic Crossover in Protein Hydration Water. Proc. Natl. Acad. Sci. USA 2006, 103, 9012–9016. (47) Pawlus, S.; Khodadi, S.; Sokolov, A. P. Conductivity in Hydrated Proteins: No Signs of the Fragile-to-Strong Crossover. Phys. Rev. Lett. 2008, 100, 108103. (48) Voronoi, G. F. Nouvelles Applications des Param`etres Continus a` la Th´eorie des Formes Quadratiques. J. Reine Angew. Math. 1908, 134, 198–287. (49) Stirnemann, G.; Laage, D. On the Origin of the Non-Arrhenius Behavior in Water Reorientation Dynamics. J. Chem. Phys. 2012, 137, 031101. (50) Fenimore, P. W.; Frauenfelder, H.; McMahon, B. H.; Young, R. D. Bulk-Solvent and Hydration-Shell Fluctuations, Similar to - and -Fluctuations in Glasses, Control Protein Motions and Functions. Proc. Natl. Acad. Sci. USA 2004, 101, 14408–14413.

20 ACS Paragon Plus Environment