How water contributes to pressure and cold denaturation of proteins

The mechanisms of cold- and pressure-denaturation of proteins are matter of debate and are commonly understood as due to water-mediated interactions. Here we study several cases of proteins, with or without a unique native state, with or without hydrophilic residues, by means of a coarse-grain protein model in explicit solvent. We show, using Monte Carlo simulations, that taking into account how water at the protein interface changes its hydrogen bond properties and its density fluctuations is enough to predict protein stability regions with elliptic shapes in the temperature-pressure plane, consistent with previous theories. Our results clearly identify the different mechanisms with which water participates to denaturation and open the perspective to develop advanced computational design tools for protein engineering.


1
Water plays an essential role in driving the folding of a protein and in stabilizing the tertiary protein structure in its native state [1,2]. Proteins can denaturate-unfolding their structure and loosing their activity-as a consequence of changes in the environmental conditions. Experimental data show that for many proteins the native folded state is stable in a limited range of temperatures T and pressures P [3][4][5][6][7][8] and that partial folding is Tmodulated also in "intrinsically disordered proteins" [9]. By hypothesizing that proteins have only two different states, folded (f ) and unfolded (u), and that the f ←→ u process is reversible at any moment, Hawley proposed a theory [10] that predicts a close stability region (SR) with an elliptic shape in the T − P plane, consistent with the experimental data [11].
High-P unfolding can be explained by the loss of internal cavities in the folded states of proteins [33], while denaturation at negative P has been experimentally observed [34] and simulated [34,35] recently. Cold-and P -unfolding can be thermodynamically justified assuming an enthalpic gain of the solvent upon denaturation process, without specifying the origin of this gain from molecular interactions [36]. Here, we propose a molecular-interactions model for proteins solvated by explicit water, based on the "many-body" water model [29,[37][38][39][40][41].
We demonstrate how the cold-and P -denaturation mechanisms can emerge as a competition between different free energy contributions coming from water, one from hydration water and another from bulk water. Moreover, we show how changes in the protein sequence affect the hydration water properties and, in turn, the stability of the protein folded state-a relevant information in protein design [25].
The many-body water model adopts a coarse-grain (CG) representation of the water coordinates by partitioning the available volume V into a fixed number N 0 of cells, each with where v 0 is the water excluded volume. Each cell accommodates at most one molecule with the average O-O distance between next neighbor (n.n.) water molecules given by r = v 1/3 . To each cell we associate a variable n i = 1 if the cell i is occupied by a water molecule and has v 0 /v > 0.5, and n i = 0 otherwise. Hence, n i is a discretized 2 density field replacing the water translational degrees of freedom. The Hamiltonian of the bulk water has a first term, summed over all the water molecules i and j at O-O distance r ij , accounting for the van der Waals interaction, with U(r) ≡ ∞ for r < r 0 ≡ v 1/3 0 = 2.9Å (water van der Waals diameter), U(r) ≡ 4ǫ[(r 0 /r) 12 − (r 0 /r) 6 ] for r ≥ r 0 with ǫ ≡ 5.8 kJ/mol and U(r) ≡ 0 for r > r c ≡ 6r 0 (cutoff).
The second term represents the directional component of the hydrogen bond (HB), with J/4ǫ = 0.3 [42], N the OOH angle is associated to a bonded state. Thus we choose q = 6 to account correctly for the entropy variation due to HB formation and breaking.
The third term, with N coop ≡ i n i (l,k) i δ σ ik ,σ il , where (l, k) i indicates each of the six different pairs of the four indices σ ij of a molecule i, accounts for the HB cooperativity due to the quantum many-body interaction [43] and leads to the low-P tetrahedral structure [44]. We choose J σ /4ǫ ≡ 0.05 ≪ J, to guarantee an asymmetry between the two HB terms.
Increasing P partially disrupts the open structure of the HB network and reduces v toward v 0 . We account for this with an average enthalpy increase P v HB /v 0 = 0.5 is the average volume increase between high-ρ ices VI and VIII and low-ρ (tetrahedral) ice Ih. Hence, the total bulk volume is HB . We assume that the HBs do not affect the n.n. distance r, consistent with experiments [44], hence do not affect the U(r) term.
Next we account for the effects of protein-water interaction. Experiments and simulations show that near a hydrophobic (Φ) residue water-water HBs are more stable then in bulk [15,[45][46][47] with stronger water-water correlation [48]. We model this by replacing J of Eq. (1) with J Φ > J for HBs at the Φ interface. This choice, according to Muller [36], ensures the water enthalpy compensation during the cold-denaturation [49].
The interaction at the Φ interface affects the hydration water density and fluctuations.
Some works suggest a decrease of interfacial water density [50][51][52][53], while recent simulations 3 show an increase of density in the first hydration shell of any solute [54] and an increase of compressibility near Φ solutes with size 0.5 nm for water [26,48,55] or water-like solvents [56] with respect to bulk. Increasing P induces a further increase of density and reduces the effect of the Φ interface on the compressibility of the hydration shell [26,48]. We incorporate this behavior into the model by using the following considerations. From the equilibrium condition for the thermodynamic potential of hydration water and the coexisting vapor at the Φ interface at fixed T , according to the Eq.
where v (Φ) is the volume per hydration water molecule and P * < 0 is the equilibrium vapor pressure at the given T . If we attribute this P -dependence to the interfacial HB properties and expand it as a power series in P , the average volume change per water-water HB at the Φ interface is where v (Φ) HB,0 is the volume change associated to the HB formation in the Φ hydration shell at P = 0, k i > 0 ∀i and lim P →∞ v where V (Φ) and N (Φ) HB are the Φ hydration shell volume and number of HBs, respectively. Because we are interested to small values of P , i.e. near the biologically relevant atmospheric pressure, we include in our calculations only the linear term in Eq.(2) [58]. In the following we fix HB /v 0 = 0.5 and J Φ /J = 1.83. Our results have minor qualitative differences by including up to the third order in Eq.(2) or by changing up to 50% the parameters.
Because our goal here is to calculate the water contribution to denaturation, we model the protein as a self-avoiding Φ homopolymer, without internal cavities [59], whose residues occupy n.n. cells with no residue-residue interaction but the excluded volume, as in other CG approaches to the problem [13][14][15]60]. This implies that the protein has several "native" states, all with the same maximum number n max of residue-residue contacts. To simplify the discussion, we initially neglect energetic contributions of the water-Φ residue interaction. 3D do not show qualitative changes, mainly because the number of n.n. water molecules is four both in 2D and 3D for the tendency of water to form tetrahedral structures in 3D.
We consider that the protein is folded if the average number of residue-residue contacts n rr ≥ 50% n max . We find an elliptic SR (Fig.1), consistent with experiments and the Hawley theory [10,11], with heat-, cold-, and P -unfolding. The elliptic shape is preserved when we change the threshold for n rr , showing that the f ←→ u is a continuous process. In the SR the folded protein (Fig. 2a) minimizes the number of hydrated Φ residues, reducing the energy cost of the interface, as expected.
First, we observe that the model reproduces the expected entropy-driven f ←→ u for increasing T at constant P (Fig. 2b). The entropy S increases both for the opening of the protein and for the larger decrease of HBs.
Upon isobaric decrease of T the internal energy dominates the system Gibbs free energy (Fig. 2c). However, N HB saturates at T lower than the SR, therefore the only way for the system to further minimize the internal energy is to increase N  Upon isothermal increase of P the protein denaturates (Fig. 2d). We find that this change is associated to a decrease of N (b) HB and a small increase of N (Φ) HB that lead to a net decrease of V at high P , as a consequence of Eqs. (2) and (3), and an increase of internal energy. At high P the P V -decrease associated to the f −→ u process at constant T dominates over the concomitant internal energy increase, determining a lower Gibbs free energy for the u state.
Hence water contribution to the high-P denaturation is density-driven, as emphasized by the increase of local density near the unfolded protein.
Finally, upon isothermal decrease of P toward negative values (Fig. 2e), the enthalpy decreases when the contribution (P v HB decreases, i.e. when N (Φ) HB increases. Therefore we find that under depressurization the denaturation process is enthalpy-driven.
From the Clapeyron relation dP/dT = ∆S/∆V applied to the SR [10], we expect that the f ←→ u process is isochoric at the SR turning points where ∂T /∂P | SR = 0, while is isoentropic at the turning points where ∂P/∂T | SR = 0. In particular, at any T and P the volume change in the f −→ u process is given by We estimate the Eq. (4) calculating the average volume V u and V f in a wide range of T and P , equilibrating water around a completely unfolded protein state and a completely folded state (with n rr = n max ). Consistently with the Hawley's theory [10], we find that the T -denaturation is accompanied by a positive entropy variation ∆S > 0 at high T and an entropic penalty ∆S < 0 at low T , while the P -denaturation by a decrease of volume ∆V < 0 at high P and an increase of volume ∆V > 0 at low P (Fig. 3). In particular, at P = 0.3(4ǫ/v 0 ), corresponding to ≈ 500 MPa, we find ∆V ≈ −2.5v 0 , hence |P ∆V | = 0.75(4ǫ) ≈ 17 kJ/mol, very close to the typical reported value of 15 kJ/mol [7]. By varying the parameters v (Φ) HB and J Φ we find that the first is relevant for the P -denaturation, as expected because it dominates Eq. (3), while the second affects the stability range in T . Both combine in a non-trivial way to regulate the low-T entropic penalty. We test our results including a small water-Φ residue attraction and find no qualitative differences but a small change in the T -range of stability of the folded protein.
Next, we study the case of a protein with hydrophobic (Φ) and hydrophilic (ζ) residues [12,61], with a residue-residue interaction matrix A i,j = ǫ rr if residues i and j are n.n.
in the unique native state, 0 otherwise. Water molecules interact with energy ǫ w,Φ < J  Fig. 1. We test that changing the parameters, within physical ranges, modifies the SR, reproducing a variety of experimental SRs [11], but preserves the elliptic shape. 8 and ǫ w,ζ > J with n.n. Φ and ζ residues respectively, accounting for the polarization of the solvent near the ζ residues. The polar ζ residues disrupt the tetrahedral order of the surrounding water molecules. Thus we assume that a ζ residue j and a n.n. water molecules i form a HB when the latter has σ i,j in the state q (ζ) j = 1, . . . , q preassigned to j. Finally, we consider that water-water enthalpy in the hydration shell is H λ,λ ≡ −J λ + P v (λ) HB , if both molecules are n.n. to the same type of residue or H λ,µ ≡ (H λ,λ + H µ,µ )/2 if the n.n. residues belong to different types, with λ, µ = Φ, ζ, and J ζ ≤ J [62] (Fig.4).
Despite the complexity of the heteropolymer model, we find results that are similar to the homopolymer case, with the qualitative difference that with ζ-residues we find a locus ∆V = 0 with negative slope and increased stability toward (i) cold-and (ii) Pdenaturation. In particular, for our specific choice of parameters, for the heteropolymer the cold denaturation at P = 0 occurs below the glass transition, instead of ≈ 50% above as for the homopolymer. Furthermore, the SR against P is ≈ 2 times larger with ζ-residues than without. This comparison suggests that the water contribution is relevant to the f ←→ u independently on the residue sequence, although the residue-residue interactions increase the stability of the folded state.
In conclusion, our model for protein folding reproduces the entire protein SR in explicit solvent and allows us to identify how water contributes to the T -and P -denaturation processes. The model is thermodynamically consistent with Hawleys theory but, in addition, allows for intermediate states for the f ←→ u process. We find that cold denaturation is energy-driven, while unfolding by pressurization and depressurization, in addition to other suggested mechanisms [33], are density-and enthalpy-driven by water, respectively. For these mechanisms is essential to take into account how the protein-water interactions affect the stability of the water-water HB and the water density in the hydration shell. In particular, both properties control the low-T entropic penalty. Our results are qualitatively robust against modification of the model parameters, within physical ranges, and the model is computationally efficient thanks to the adoption of a CG water model, representing a step towards the development of a theoretical and computational approach for protein design and engineering.