After decades of using urea as denaturant, the kinetic role of this molecule in the unfolding process is still undefined: does urea actively induce protein unfolding or passively stabilize the unfolded state? By analyzing a set of 30 proteins (representative of all native folds) through extensive molecular dynamics simulations in denaturant (using a range of force-fields), we derived robust rules for urea unfolding that are valid at the proteome level. Irrespective of the protein fold, presence or absence of disulphide bridges, and secondary structure composition, urea concentrates in the first solvation shell of quasi-native proteins, but with a density lower than that of the fully unfolded state. The presence of urea does not alter the spontaneous vibration pattern of proteins. In fact, it reduces the magnitude of such vibrations, leading to a counterintuitive slow down of the atomic-motions that opposes unfolding. Urea stickiness and slow diffusion is, however, crucial for unfolding. Long residence urea molecules placed around the hydrophobic core are crucial to stabilize partially open structures generated by thermal fluctuations. Our simulations indicate that although urea does not favor the formation of partially open microstates, it is not a mere spectator of unfolding that simply displaces to the right of the folded←→unfolded equilibrium. On the contrary, urea actively favors unfolding: it selects and stabilizes partially unfolded microstates, slowly driving the protein conformational ensemble far from the native one and also from the conformations sampled during thermal unfolding.
Urea is a protein denaturant that has been used for decades in the study of protein folding/unfolding; however, after many years of research the ultimate reasons of the denaturing properties of urea remain elusive [1], [2]. The dominant paradigm for unfolding (the “direct” mechanism) claims that the denaturant properties of urea are related to its capacity to interact with exposed protein residues more strongly than water [3]–[15]. However, the nature of such a preferential interaction is not so clear. Thus, while some authors suggest that it is mostly electrostatic and related to the formation of direct hydrogen bonds [7]–[9], [16]–[17], others claim that preferential dispersion is the leading term [13]–[15]. It is also unclear whether the major destabilizing effect of urea is related to interaction with the backbone [6]–[7] or with side chains [8]–[12]. In the latter case, there is also discussion regarding the preferential side chains: polar and charged [9] or apolar [4], [10]–[12]. We recently combined multi-replica molecular dynamics (MD) simulations and direct NMR measures of ubiquitin to characterize the “urea unfolded ensemble” of this model protein [15]. Our results suggest that urea stabilizes flexible over-extended conformations of the protein, which are unlikely to be sampled in the “unfolded” state of aqueous proteins. Extended conformations of the protein with exposed hydrophobic surfaces are more urea-philic than the native globular state, due mostly to extensive London dispersion interactions (the attractive contribution in Van der Waals interactions between instantaneous dipoles) between apolar side chains and urea molecules in the first solvation shell of unfolded conformations. We believe that our results in reference 15 clarify the molecular basis of the effect of urea on the thermodynamics of the folded←→unfolded equilibrium, but unfortunately, they do not provide information on the kinetic role of urea in the unfolding process. In other words: does urea actively induce protein unfolding? Or, on the contrary, does it passively stabilize the unfolded state by selectively binding to unfolded conformations? To analyze this point, we should characterize the effect of urea in the first stages of thermochemical unfolding, when the protein structure is still close to the native conformation and internal residues are not fully exposed. Clearly, a study of this nature presents many difficulties, the most important being that the effect of urea on early stages of unfolding might be dependent on the native structure. Therefore, to obtain conclusions of general validity, all representative protein folds should be addressed. Also, results can be force-field-dependent, so if we aim to obtain robust conclusions, we should perform simulations with a variety of force-fields. Given the typical kinetics of the folding/unfolding transitions of small globular proteins [18], microsecond (µsec) long simulations should trace the first stages of these processes. In the current work, we investigate the first stages of urea-driven protein unfolding using µsec-long atomistic simulation; to gain universality, we used 30 proteins representative of all protein folds, while to protect our conclusions from force-field-related uncertainties, we used several of the most popular force-fields. The results derived from this study provide a robust and complete picture of the role of urea in destabilizing folded states of proteins, and more importantly, on the molecular mechanisms by means of which urea contributes to accelerating protein unfolding.
We first validated our protocol using three ultra-representative proteins (in bold in Table 1), one for each of the main classes in the Structural Classification of Proteins (SCOP, [19]). We monitored the protein stability in three environments: i) in chemical unfolding conditions, in 8M urea and with a mildly high temperature (T = 368K) to speed up the observable effects; ii) in thermal unfolding condition, in water with the same high temperature; this control allowed us to distinguish the effect of urea and temperature on protein unfolding; iii) in water at room temperature as final control. Four force-fields were used (OPLSAA - ON2; CHARMM - C22; AMBER99 - P99 and P99SBILDN) for each system (see Methods for the description of the force-fields used), collecting in total 36 simulations of 1-µsec length each. After the validation of our protocol, we extended the chemical unfolding simulations to a larger set of proteins, to avoid any bias in the conclusion due to the native structure. We performed 1 µsec of simulation in urea at high temperature (T = 398K) for 30 proteins covering all the major protein folds (Table 1 and Suppl. Dataset S1). Each system was simulated in three force-fields (C22, ON and P99), excluding P99SBILDN as reported above, and collecting in total 90 simulations. To have a more realistic picture of the native state, instead of using the crystal structure, we used as control 0. 1 µsec-long simulations in water at room temperature for all the 90 systems. The analysis described here reveals some common robust trends that illustrate the effect of urea during the early stages of protein unfolding.
MD simulations with additive potentials and explicit solvent have become very popular to explore chemical unfolding of protein. There is little doubt that the use of the technique has produced sizeable advances in the field, but we cannot ignore some potential caveats in the beginning of this discussion. First, for computational reasons we (and most authors in the field) are using classical non-polarizable force-fields, which might not be accurate enough to deal with a complex process such as unfolding. Previous studies [4]–[5], [7]–[15], [25] have however demonstrated that urea/water/protein effective parameters are able to reproduce a variety of experimental observables, such as mass densities and radial distribution functions of urea/water solutions derived from neutron scattering experiments [25], the experimental water/urea transfer free energies of tripeptides [10], and the urea density around unfolded proteins found by vapor pressure osmometry measures [4], [8], [10], [13], [15], [24]. Furthermore, our recent work [15] has demonstrated that unbiased MD simulations in 8M urea reproduce very accurately the unfolded ensemble as determined from a variety of spectroscopic techniques (including SAXS and NMR) under the same conditions. Thus, despite their simplicity current force-fields reproduce reasonably well urea/water/protein mixtures. We should remember that since we are exploring a microsecond-long process, no direct experimental data is available for comparison and accordingly caution is required. This move us to use a consensus approach, running the simulations with different force-fields to extract those results that seem robust to force-field changes. A second reason of concern is related to the stochastic nature of unfolding, where individual trajectories can show different degree of unfolding [13]. Again, by comparing different trajectories we tried to define robust findings, but we cannot ignore that the experimental result is the averaging a near-Avogadro number of trajectories. A third reason of concern, common to many experimental studies, is the generality of the results, i. e. how general are the results obtained with a few model proteins. To convince ourselves on the general validity of our results we repeated the unfolding studies for a large number of proteins representative of all prevalent folds. Despite the obvious caveat of any theoretical study, this approach provided a picture of unprecedented, to our knowledge, completeness and robustness of the early stages of urea unfolding. Under our simulation conditions (8M urea at T = 368 K), we detected clear signals of unfolding in the microsecond range, but progress in denaturation was smaller than that reported for model proteins of reduced stability [13]. Overall, the advance in unfolding of proteins at T = 368 after 1 µsec of MD is not dependent on the fold, nor secondary structure composition, and was similar for proteins with and without disulphide bridges in the native form. No dramatic differences in the advance of unfolding were found between water and urea simulations performed at the same temperature; however, unfolding paths in the presence or absence of urea differed, since the partially unfolded structures sampled in hot water maintained the hydrophobic residues hidden in the core of the protein, while such residues were more accessible in presence of urea. Urea solutions are more viscous than pure water, which in our simulation reduced high-frequency movements in the protein, generating an unexpected slow down of the atomic-motions. Urea residence times around protein residues were large, especially when urea molecules diffuse close to the hydrophobic core or to the interface between rigid and thermally mobile regions (hinge points). We consider that the sticky nature of urea and its preferential placement at hinge points is crucial for unfolding, since it favors the rapid trapping of residues that become exposed as a consequence of stochastic thermal motions. The stabilizing effect of urea on exposed residues slowly biases the trajectory towards the unfolded state, by decreasing the chances of microscopic refolding [12], [26]. The effect of stabilization of exposed residues is especially productive in terms of unfolding when residues are apolar, since in this case urea (but not water) traps very efficiently the residue, increasing the accessibility of other apolar residues in the vicinity. The ensuing greater recruitment of urea in the region leads to a cooperative effect resulting in the acceleration of protein unfolding. Our data shows that, similar to the unfolded state [13], [15], it is the van der Waals interactions that drive the accumulation of urea on the surface of the folded protein. However, the role of H-bonding cannot be dismissed, as these bonds are crucial for the stabilization of long-living urea interactions near hinge points, which in turn are required to bias intrinsic protein dynamics towards unfolding. Clearly, “direct” effects not only are the main factors responsible for the urea-mediated stabilization of the unfolded state [15], but are also relevant in guiding the first steps of urea unfolding. Microscopic unfolding events are related to stochastic thermal motions, which are in principle similar to those that occur spontaneously in water at room temperature. However, urea is not a mere passive spectator that simply stabilizes the small percentage of unfolded protein coexisting within the native ensemble and leading to a displacement in the folded←→unfolded equilibrium towards the denatured state. On the contrary, urea has a dual function: i) it takes advantage of microscopic unfolding events, decreasing their chances of refolding, and favoring further unfolding [12], [26]; and ii) among these microscopic unfolding events it selects and stabilizes microstates with exposed hydrophobic regions [4] (see Suppl. Figure S7). These effects lead to a slow divergence in the temperature-unfolding pathways in water and urea, and, as shown for ubiquitin [15], to distinct unfolded states. Consequently, concepts such as folded and unfolded states or folding and unfolding pathways need to be revisited and reformulated considering the nature of the denaturant used.
As model structures for the main protein-folds we used the same structures selected in our previous work in reference 20. We first explored the early stages of urea unfolding using three ultra-representative proteins for the most populated fold in the three main classes in the SCOP database (all-α 1OPC, all-β 1CQY and α/β 1KTE; [19], [20]). Once the simulation protocols had been validated with these proteins, the study was extended to a larger set, consisting of 30 structures (110 residues on average) representative of the most populated protein folds ([20], [27] and Suppl. Dataset S1) All starting structures were taken from the Protein Data Bank (PDB; [28]) and processed using our standard procedure implemented in the MDWeb server [29]: experimental structures were titrated to define the major ionic state at neutral pH, neutralized by ions (sodium and chloride), minimized for 1000 steps, heated up to the final temperature, and solvated using a 8M urea/water octahedron box with a spacing distance of 15 Å around the system. The box was previously equilibrated in a Monte Carlo simulation using the BOSS program [30]. The water model was taken from Jorgensen' s TIP3P [31], while ion and urea force-field parameters were those considered as the default of each force-field. Urea parameters from Smith et al. [32] were used for OPLS and P99SBILDN simulations, the same charges but scaled according to the amber force-field were used in PARM 99, while Nilsson' s parameters were used in the CHARMM 22 force field [33]. Systems were then pre-equilibrated for 0. 5 ns with parm99-AMBER force field in keeping the backbone restrained by intra-molecular harmonic potentials and then equilibrated (0. 5 ns) in each force field parameters removing backbone constraints. For the small set of ultra-representative proteins, three sets of simulations corresponding to water at room temperature (T = 300 K), hot water (T = 368 K), and urea at high temperature (T = 368 K) were carried out. For each condition, we performed 1 µsec simulations using four force-fields: three general purpose ones (OPLSAA -ON2- [34]; CHARMM -C22- [35]; AMBER99 -P99- [36]), and a last-generation force-field able to accurately reproduce folded proteins (P99SBILDN, [37]). For the extended set of 30 proteins, control simulations in water were limited to 0. 1 µsec at room temperature, while the 8M urea simulations were performed, as above, for 1 µsec at T = 368 K. Simulations for the extended set of proteins were carried out using ON2, C22 and P99. All simulations were performed using periodic boundary conditions and particle Mesh Ewald [38] corrections for the representation of long-range electrostatic effects using a 1. 0 Å grid spacing and a 9 Å cutoff. All trajectories were collected with the NAMD2 [39] program. Integration of equations of motions was performed every 2 fs after removing vibrations of bonds involving hydrogen atoms using SHAKE/RATTLE algorithm [40], [41]. All simulations were carried out in the isothermal (T = 300 or 368 K) /isobaric ensemble (P = 1 atm) using the Langevin thermostat and barostats [42], [43]. The trajectories were analyzed using VMD [44] and the MdWeb server [29], as well as Flexserver which can be accessed at: http: //mmb. pcb. ub. es/FlexServ/ (see also Suppl. Text S1 for a detailed explanation of the metrics).
The delicate equilibrium between the folded and functional structure of a protein and its unfolded state is highly dependent on environmental variables such as the solvent. For example the co-solvent urea is a well-known protein denaturant that displaces the equilibrium towards unstructured and non-functional conformations of proteins. However the molecular mechanism behind its ability remains an enigma and the interpretation of the experimental data is still ambiguous. By analyzing a set of representative proteins through extensive molecular dynamics simulations in urea, we provide a robust and consensus picture of the first stages of urea-driven protein unfolding and elucidate the role of urea in accelerating protein unfolding. Our results suggest that urea, thanks to its stickiness and slow diffusion, benefits from the intrinsic flexibility of proteins and stabilizes partially open-states, slowly driving the protein toward unfolding.
Exploring Early Stages of the Chemical Unfolding of Proteins at the Proteome Scale