Role of Electrostatic Interactions in Binding of Peptides

Role of Electrostatic Interactions in Binding of Peptides and Intrinsically Disordered Proteins to Their Folded Targets: 2. The Model of Encounter Complex Involving the Double Mutant of the c‑Crk N‑SH3 Domain and Peptide Sos ABSTRACT: In the first part of this work (paper 1, Xue, Y. et al. Biochemistry 2014, 53, 6473), we have studied the complex between the 10-residue peptide Sos and N-terminal SH3 domain from adaptor protein c-Crk. In the second part (this paper), we designed the double mutant of the c-Crk N-SH3 domain, W169F/Y186L, with the intention to eliminate the interactions responsible for tight peptide−protein binding, while retaining the interactions that create the initial electrostatic encounter complex. The resulting system was characterized experimentally by measuring the backbone and side-chain 15N relaxation rates, as well as binding shifts and 1 N H temperature coefficients. In addition, it was also modeled via a series of ∼5 μs molecular dynamics (MD) simulations recorded in a large water box under an Amber ff99SB*-ILDN force field. Similar to paper 1, we have found that the strength of arginine-aspartate and arginine-glutamate salt bridges is overestimated in the original force field. To address this problem we have applied the empirical force-field correction described in paper 1. Specifically, the Lennard-Jones equilibrium distance for the nitrogen−oxygen pair across Arg-to-Asp/Glu salt bridges has been increased by 3%. This modification led to MD models in good agreement with the experimental data. The emerging picture is that of a fuzzy complex, where the peptide “dances” over the surface of the protein, making transient contacts via salt-bridge interactions. Every once in a while the peptide assumes a certain more stable binding pose, assisted by a number of adventitious polar and nonpolar contacts. On the other hand, occasionally Sos flies off the protein surface; it is then guided by electrostatic steering to quickly reconnect with the protein. The dynamic interaction between Sos and the double mutant of c-Crk N-SH3 gives rise to only small binding shifts. The peptide retains a high degree of conformational mobility, although it is appreciably slowed down due to its (loose) association with the protein. Note that spin relaxation data are indispensable in determining the dynamic status of the peptide. Such data can be properly modeled only on a basis of bona fide MD simulations, as shown in our study. We anticipate that in future the field will move away from the ensemble view of protein disorder and toward more sophisticated MD models. This will require further optimization of force fields, aimed specifically at disordered systems. Efforts in this direction have been recently initiated by several research groups; the empirical salt-bridge correction proposed in our work falls in the same category. MD models obtained with the help of suitably refined force fields and guided by experimental NMR data will provide a powerful insight into an intricate world of disordered biomolecules. I n very broad terms, proteins and nucleic acids that function as a part of cell machinery should be characterized from two different perspectives: structural and dynamic. In the field of structural studies, the preeminent technique is X-ray crystallography. In the area of dynamics, the best results are obtained through a combination of experimental NMR spectroscopy and molecular dynamics (MD) simulations. In particular, this is true for partially or fully disordered proteins, which comprise a significant and functionally important part of the eukaryotic proteome. In this area, MD simulations have been used in conjunction with the NMR data to generate conformational ensembles representative of disordered protein states.1−6 Furthermore, in the increasing number of cases, the MD trajectories are viewed as bona f ide models of disordered Figure 1. (A) Peptide-binding interface of the c-Crk N-SH3 domain as seen in the crystal structure 1CKB.16 Protein atoms within 3.9 Å of the hydrophobic (PPPVPP) portion of the bound Sos peptide are colored green. Carboxyl oxygen atoms from contiguous residues D147, E148, E149, D150, and E166, E167 are colored red. Note that there is a significant amount of local mobility associated with solvent-exposed Glu and Asp side chains.23 (B) Schematic free energy landscape of the binary system consisting of Sos and wtSH3 domain (blue line) or, alternatively, Sos and dmSH3 (orange line). The electrostatic encounter complex Sos·SH3 is comprised of many conformationally diverse states. In the plot it is separated from the stable complex Sos:SH3 by a certain low-amplitude barrier. Indeed, it appears that the salt-bridge interactions observed in Sos·wtSH3 may not be fully compatible with the hydrophobic packing and hydrogen bonding characteristic of Sos:wtSH3, i.e., the two states are separated by a barrier.12 proteins (peptides) and used to calculate NMR observables that depend on temporal variables (i.e., motional correlation times).7−12 This paper is the second part of the two-part report where we seek to characterize the electrostatic encounter complex between the disordered peptide and its folded target protein. In the first part, subsequently referred to as paper 1,12 we have examined the binding of 10-residue peptide from Ras activator protein Sos to the N-terminal SH3 domain from adaptor protein c-Crk. We were able to reproduce the electrostatically driven binding of the peptide to the SH3 domain in microsecond-scale explicit-solvent MD simulations using Amber ff99SB*-ILDN force field.13−15 The simulation starts with the free peptide (Sos) and protein (wtSH3) separated by a distance of more than 10 Å. Under the pull of long-range electrostatic forces, the peptide descends to the surface of the protein and forms the so-called electrostatic encounter complex (Sos·wtSH3), a highly dynamic state where Sos is loosely anchored to SH3 via a number of transient salt bridges. It is then a matter of time for the peptide to find a correct binding pose and form a stable complex (Sos:wtSH3). In our simulations we were able to reproduce the coordinates of the bound peptide to within ∼2 Å of the crystallographic structure.16 The binding parameters of the Sos:wtSH3 complex under the experimental conditions chosen in our study are Kd = 1.5 μM, kON = 1.5 × 109 M−1 s−1. This means that the complex has a sufficiently long lifetime, ∼0.5 ms. Indeed, once formed, Sos:wtSH3 does not dissociate during the course of several microsecond long MD trajectories. This also means that one can easily prepare an NMR sample where the protein is saturated with the peptide or vice versa. For example, in our experimental work, we used a sample where 99.9% of the (isotopically labeled) peptide was in the bound form. The stable nature of Sos:wtSH3 makes it possible to directly characterize this complex via experimental measurements and computer simulations. The comparison of the experimental NMR results with the predictions from MD simulations exposed one problem area: as it turns out, the simulations tend to exaggerate the strength of Arg-to-Asp/Glu salt bridges between the protein and the peptide. To further investigate this problem, we have conducted additional MD simulations using several popular force fields, as well as different water models, etc. As it appears, none of them can accurately model these particular salt bridge interactions, the observation that is also confirmed by recent small-molecule studies.17,18 To remedy this problem we have introduced an empirical correction into the Amber ff99SB*-ILDN force field. Specifically, we have increased by 3% the Lennard-Jones equilibrium distance for the pairwise interaction between the nitrogen atoms in Arg guanidine group and oxygen atoms in Asp (Glu) carboxyl group. This correction, termed λ = 1.03, was a result of ad hoc optimization against the experimental data, as detailed in paper 1. Note that this modification has a highly selective character and does not influence any of the Arg or Asp (Glu) interactions except in a situation when the two respective side chains come into direct contact leading to a formation of the salt bridge. Subsequent microsecond timescale MD simulations using the modified force field confirmed good agreement with the experimental NMR data, including 15 N backbone and side-chain spin relaxation data.12 These results are consistent with recent observations of Best and coworkers, who concluded that current force fields tend to overestimate the strength of protein−protein interactions, in particular, that of transient weakly specific interactions.19,20 Thus, in paper 1 we have obtained a self-consistent description of the Sos:wtSH3 complex, including transient interactions between the arginine-rich tail of Sos and the acidic patch on the surface of wtSH3. Furthermore, in our MD simulations we have observed a short-lived binding intermediate corresponding to the electrostatic encounter complex Sos·wtSH3. However, experimental characterization of Sos· wtSH3 remains beyond reach. Indeed, this state is too sparsely populated and short-lived to allow for direct characterization by NMR experiments. Although a number of NMR techniques are sensitive to low-population species (e.g., relaxation dispersion measurements, paramagnetic relaxation enhancement measurements, etc.), they are generally effective only when the population of minor species exceed ∼1%. One should also bear in mind that these experiments are indirect, i.e., the Figure 2. In silico scan mutagenesis of c-Crk N-SH3 domain: Y186a and W169a/Y186L, where a is a variable amino acid (left and right portions of the plot, respectively). Shown are the predicted changes in the binding free energies for the complex of Sos with SH3 domain and the changes in thermodynamic stability of SH3 (blue and red bars, respectively). The data are from the program FoldX,24,25 which has been used with default parameters. information about minor species is gleaned from the data that are collected from the dominant species. Translating such experimental data into structural models representative of the minor states is a challenging problem.21 The level of difficulty is further increased when a three-state exchange scheme needs to be taken into consideration, Sos + wtSH3 ↔ Sos·wtSH3 ↔ Sos:wtSH3. Under such circumstances, recent relaxation dispersion study of the Fyn SH3 domain binding a 12-residue peptide failed to produce any direct evidence of the encounter complex.22 In this situation we resort to an alternative strategy to characterize the electrostatic encounter complex between Sos and SH3. Shown in Figure 1A is the peptide-binding interface of the c-Crk N-SH3 mapped according to the crystallographic coordinates of the Sos:wtSH3 complex.16 The contact area associated with the hydrophobic N-terminal portion of Sos is painted green. This region is responsible for the (relatively) tight peptide binding. The Glu/Asp carboxyl groups forming a negatively charged patch on the surface of the protein are colored red. This region interacts with the arginine-rich Cterminal tail of Sos. It is responsible for the uncommonly high on-rate kon on the path to the electrostatic encounter complex. We have designed two SH3 mutations, Y186L and W169F, that attenuate (abrogate) tight peptide binding, while preserving the electrostatic (Asp/Glu)-to-Arg interactions. The idea is to render Sos:SH3 thermodynamically unfavorable so as to shift the balance toward the Sos·SH3 state. This is illustrated schematically in Figure 1B. The blue line in this figure illustrates the energy landscape of the original system, where the global minimum corresponds to the stable Sos:SH3 complex. The orange line represents the energy profile changed by the two point mutations in the SH3 domain. The mutations interfere with the binding of the hydrophobic N-terminal region of the peptide (cf. Figure 1A), so that the high-affinity complex becomes disrupted. On the other hand, the mutations do not directly affect the Asp/Glu-rich patch. Hence the part of the phase space corresponding to the electrostatic encounter complex Sos·SH3 remains unchanged. Figure 1B illustrates the conceptual approach to mutant design. The picture is clearly idealizedwe do not expect the complex between Sos and the double mutant SH3 (dmSH3) to be a perfect representation of the original electrostatic encounter complex. Nevertheless we suggest that it provides a reasonably good model for the said encounter complex. Hence we refer to it as Sos·dmSH3. Importantly, this engineered state is strongly populated, which allows us to directly characterize it via 15N spin relaxation rates and other conventional experimental parameters. ■ MATERIALS AND METHODS Mutant Design. The scanning mutagenesis of c-Crk N-SH3 was first conducted in silico using the program FoldX.24,25 Briefly, we evaluated the effect of site mutations at the hydrophobic interface (colored green in Figure 1A) searching for mutations that disfavor peptide binding but do not affect the stability of the protein. For example, the mutation Y186L is predicted to reduce the binding affinity of Sos by 1.9 kcal/mol, while at the same time stabilizing protein structure by −0.4 kcal/mol (cf. dashed line in the left half of Figure 2). Furthermore, double mutation W169F/Y186L reduces the binding affinity of Sos by 3.2 kcal/mol relative to the wild-type protein, while stabilizing protein structure by −0.5 kcal/mol Figure 3. (A) 1HN-15N HSQC titration of 15N-labeled Sos with unlabeled dmSH3. The starting point corresponds to 0.24 mM Sos, no dmSH3 (red contours in the plot); the final point corresponds to 0.13 mM Sos, after accounting for dilution, and 1.22 mM dmSH3 (blue contours in the plot). The minor peaks, labeled with asterisks, represent distinct conformational species arising from cis−trans isomerization of prolines. Unlike in the case of wtSH3 binding, these minor peaks move during the course of the titration. This can be understood considering that the PPPVPP portion of the peptide is no longer the determinant of binding and, therefore, the binding affinity of the minor species, Sos*, becomes comparable to that of the dominant species, Sos. (B) Concentration dependence of chemical shifts from 1HN-15N HSQC titration of Sos with dmSH3. The data have been fitted on per-residue basis using fast-exchange model, δ = pfreeδfree + pboundδbound, where pfree is the fraction of the free Sos, pfree = ((cSos − cSH3 − Kd) + ((cSos − cSH3 − Kd)2 + 4cSosKd)1/2)/(2cSos), pbound is the fraction of the bound Sos, pbound = 1 − pfree, δfree and δbound are respective chemical shifts, cSos and cSH3 are the analytical concentration of the protein and the peptide, and Kd is the dissociation constant of the complex. The fitting procedure involves two fitting parameters per residue: Kd and δbound. The extracted values of Kd are 192, 213, 167, 323, and 234 μM for R8, V4, R7, R9, and R10, respectively (listed in the order of the titration curves in the plot, from the top down). The titration data under fast-exchange conditions do not allow for extraction of kON rates (cf. paper 112). (dashed line in the right half of Figure 2). While these estimates are, of course, approximate, they allow us to approach the process of mutant design in a rational manner. Of note, using the Kd value of Sos:wtSH3 (see above) we can easily estimate the expected binding affinity of W169F/Y186L mutant, ∼325 μM. As it turns out, this is very close to the experimentally determined value (see below). We assume that the drop in the binding affinity is mainly due to the loss of two prominent protein-to-peptide hydrogen bonds, Y186-P2 and W169-P5, that are found in the crystal structure and persist in the MD simulation of Sos:wtSH3. Experimentally, mutagenesis was conducted using a QuikChange 2.0 kit. The 1HN-15N HSQC spectrum of W169F/Y186L SH3 (dmSH3) superimposed onto the spectrum of wtSH3 is shown in Figure S1. The spectral assignment has been verified and augmented by means of the HNCACB experiment.26 The spectral data confirm that the structure of the protein remains intact upon the two mutations. In addition to dmSH3, we have also prepared several other mutants of c-Crk N-SH3 and titrated the corresponding 15Nlabeled samples with Sos. In doing so we have (crudely) estimated that the binding affinity of Y186L and W169F/ Y186W amounts to ∼100 μM, whereas the binding affinity of F141W/W169F/Y186L reaches ∼200 μM. None of these mutants have been investigated any further, as we have focused on the dmSH3. Experimental Procedures. In addition to dmSH3 and Sos· dmSH3, we have also collected 15N relaxation data from the sample of free Sos peptide. The concentration of this sample was 0.5 mM, in line with the other samples containing labeled Sos. In the spectra of free Sos, the four arginine side-chain resonances 1Hε,15Nε were all overlapped, giving rise to a single spectral peak. Therefore, the effective average relaxation parameters have been obtained for these sites. In the spectra of Sos·dmSH3, the same resonances gave rise to two poorly resolved peaks with similar intensities. Since we could not assign those (near-degenerate) peaks, we also chose to determine the effective average relaxation rates. In doing so we fitted the spectral data assuming that there are two partially overlapped peaks, then added the two obtained peak intensities and used the results to determine average relaxation parameters. Temperature coefficients of backbone 1 H N chemical shifts27 were obtained from HSQC spectra recorded at 5 and 23.4 °C. All other experimental details are the same as described in paper 1. Kd for Binding of Sos to dmSH3. NMR titration experiment has been carried out by adding unlabeled dmSH3 to the sample of 15N-labeled Sos, see Figure 3A. Compared to the wild-type protein, the titration effects are far less pronounced. Nevertheless, the results suggest that the peptide and the protein engage in a meaningful interaction. In part, the shifts are reduced because the ring-current effects across the binding interface are decreased as a consequence of the mutations in the protein, Y186L and W169F. More significantly, though, Sos·dmSH3 represents a “fuzzy” complex where Sos retains a significant motional freedom while coordinated to the surface of dmSH3 (see below). This leads to partial averaging of the shifts such as the 1HN shift in residue R8. The titration data from all Sos residues are consistent with fast exchange between the two inequivalent sites (free and bound forms of Sos). The corresponding fits are shown in Figure 3B leading to the consensus estimate Kd ∼ 200 μM. It is important to keep in mind that “bound Sos” in fact represents the ensemble of rapidly interconverting conformational species. Tumbling Time of dmSH3 and Sos·dmSH3. Using 15N R2/R1 data, we have determined the tumbling time of dmSH3 at 2 mM protein concentration. The result, τrot = 5.8 ns, is appreciably lower than for wtSH3 under the same conditions, yet significantly higher than in the dilute sample (see Table S1). The anisotropy and the orientation of the long axis of the rotational diffusion tensor also experience a noticeable change compared to the wtSH3. All of this leads us to conclude that the two mutations, Y186L and W169F, substantially reduce the self-association effect in c-Crk N-SH3 without completely eliminating it (see paper 112). This finding is not surprising. Both Y186 and W169 are among the 10 residues that show the strongest chemical shift perturbation upon dimerization. At the same time, there are other dimerization sites that are widely distributed over the surface of the protein. As a next step, we have prepared the sample containing 2 mM dmSH3 and 0.5 mM Sos (these are the sample conditions subsequently used for dynamics studies). Under these conditions, 22% of the protein is loaded with the peptide, increasing the apparent tumbling time of the protein from 5.8 to 6.1 ns. This result can be readily extrapolated to estimate the correlation time of Sos·dmSH3 complex, τrot = 7.1 ns. Unlike in the case of Sos:wtSH3, we cannot rely on the peptide data to verify this result (as discussed later, Sos retains a significant amount of motional freedom when coordinated to dmSH3). Therefore, we prepared an additional control sample containing 2 mM dmSH3 and 7.5 mM Sos. Under these (nearsaturating) conditions, 97% of the protein is loaded with the peptide. The tumbling time τrot measured in this sample turns out to be 6.7 ns. After correcting for a small fraction of the free protein, we arrive at the value 6.8 ns, which is in good agreement with the previous estimate, 7.1 ns. Finally, we assume that the tumbling time of Sos·dmSH3 complex in our relaxation experiments is 7.0 ns, which should be accurate to within ∼3%. This is the same as previously determined for Sos:wtSH3 (see paper 112). Correcting Sos·dmSH3 Data to Account for Free Sos. Using the dissociation constant determined for Sos·dmSH3 complex, Kd = 200 μM, we can readily estimate that in our working sample (0.5 mM Sos, 2.0 mM dmSH3) only 89% of the peptide is bound to the protein while the remaining 11% are free. The proportion of the free peptide is non-negligible and should be taken into consideration when interpreting the results of the relaxation measurements. To begin with, it can be readily demonstrated that the exchange between Sos·dmSH3 and free Sos is fast on the NMR time scale. Recall that the onrate for binding of Sos to c-Crk N-SH3 is determined by strong electrostatic attraction. The mutations in dmSH3 are designed such as to not interfere with electrostatic interactions. Therefore, it is reasonable to suggest that the kON rate in Sos·dmSH3 sample remains very similar to the one found in Sos:wtSH3, kON = 1.5 × 109 M−1 s−1. Using this assumption, we can estimate the off-rate for Sos·dmSH3, kOFF = KdkON = 3 × 105 s−1, and further calculate the correlation time for exchange between free and protein-bound Sos, τex = (kON[dmSH3]free + kOFF)−1 ∼ 0.4 μs. This estimated τex value corresponds to the extreme narrowing regime, where the signals from Sos should experience no line broadening. Indeed, the titration spectrum in Figure 3A shows no evidence of exchange broadening. This result can be further validated by comparing the experimental 15 N R2 rates with the dipolar-CSA cross-correlated relaxation rates ηxy. While R2 rates are generally sensitive to the microsecond time-scale exchange, the cross-correlations are free from the chemical exchange effects.28 The inspection of the experimental data in Figure 6 demonstrates that R2 and ηxy profiles are similar and thereby confirms the absence of any appreciable Rex contribution. Under the fast exchange conditions, the experimentally measured Sos relaxation rates correspond to a simple weighted average, e.g., R1 = pfreeRfree + pboundRbound . Since our focus is 1 1 primarily on the bound form, we have postprocessed the experimental data to extract Rbound . Specifically, we have used 1 the experimentally measured Rfree 1 rates from the sample of free Sos and the experimentally determined populations pbound = 0.89, pfree = 0.11 to arrive at Rbound . The same method was 1 employed for all other relaxation parameters, i.e., R2, ηxy, and NOE, as well as the peptide binding shifts. Note that the results show very little sensitivity to small variations in Kd (see above). To test the validity of this correction procedure, we have prepared a control sample containing a much larger amount of the peptide (3.6 mM Sos, 2.0 mM dmSH3). In this sample, only 50% of Sos is expected to be complexed with dmSH3. After application of the same processing scheme as described above, we have obtained a new set of Rbound and Rbound rates. 1 2 These new data proved to be within 3% of the previously determined values, thus demonstrating that we can successfully remove the contribution due to free Sos. MD Simulations. The simulations have been conducted using the Amber ff99SB*-ILDN force field13−15 and TIP3P solvent29 as described in paper 1.12 Both Amber 11 and 14 packages have been employed.30,31 The potential associated with dihedral angle ψ has been corrected according to Best and Hummer;14 the correction has been applied to all residues, including Gly and Pro, same as in paper 112 (we have verified that the outcome does not depend on whether the correction is extended to glycines and prolines). Under Amber 11, LennardJones equilibrium distances for Nη(Arg), Oε(Glu) and Nη(Arg), Oδ(Asp) pairs were scaled as described in paper 1.12 Under Amber 14, the corresponding changes to the force field were made in a simpler fashion using ParmEd program within the AmberTools suite. Specifically, the original equilibrium distance r0 = 3.4852 Å was multiplied by the scaling factor λ = 1.03 using the following ParmEd input script: Listed in Table 1 are the parameters of Sos·dmSH3 trajectories discussed in this study. In addition, we have also used the 2.1 μs trajectory of free Sos, which serves as a reference to calculate binding shifts. NMR observables have been computed on the basis of the MD trajectories as described in paper 1.12 In calculating the dipolar correlation functions g(τ) for free Sos peptide, we did not superimpose the peptide coordinates and did not multiply the extracted g(τ) functions by exp(−τ/τrot). Instead we rescaled the time axis of g(τ) by a factor of 2.7 to correct for the intrinsically low viscosity of the TIP3P water model.33,34 This approach is generally applicable to disordered proteins and peptides where it is impossible to separate the overall motion from internal conformational dynamics. We have also tested this strategy for the Sos peptide in a fuzzy complex with dmSH3. Finally, the content of the helical motif in Sos was quantified using the angle ψ in residue P6. Table 1. Parameters of the Sos·dmSH3 Trajectories nomenclaturea trajectory length (μs) protein−protein separation (Å)b MD(xray,λ = 1.00) MD(xray,λ = 1.03) MD(rand,λ = 1.03)c 5.4 5.2 6.0 38.3 38.5 40.6 a The nomenclature refers to the starting coordinates of the MD trajectory and the scaling factor λ applied to the Lennard-Jones equilibrium distances for the Nη, Oε/δ atom pairs. The xray designation indicates that the starting conformations is modeled after the crystallographic structure of Sos:wtSH3 complex (1CKB), where two point mutations, Y186L/W169F, have been introduced using the LEaP facility of Amber. The rand designation indicates that the initial conformation of the peptide and its placement are random (see paper 112 for details). bMinimal interatomic distance between the simulated dmSH3 molecule and its periodic image (i.e., minimal separation). The values are obtained from the starting MD coordinates (first frame) using the TopoTools module in VMD.32 The xray coordinates were solvated with a 14 Å layer of TIP3P water; the rand coordinates, where the peptide is positioned at a certain distance from the surface of the protein, were solvated with a 10 Å layer. In our previous simulations of the more compact Sos:wtSH3 complex, the respective settings were 8 and 10 Å (these values have been indicated incorrectly in paper 112). cAn additional 6.5 μs trajectory has been recorded beginning with a random configuration; the minimum protein−protein separation in the starting frame of this additional trajectory was 35.8 Å. ■ RESULTS AND DISCUSSION Simulation of Sos·dmSH3 Using Amber ff99SB*-ILDN: Comparison with Experimental Data. In paper 1,12 we have described MD simulations of Sos peptide in a moderately high-affinity (∼1 μM) complex with the wild-type N-SH3 domain from the adapter protein c-Crk. In the simulations, we have observed that the arginine-rich C-terminal tail of Sos made exceedingly stable salt bridges with the Glu/Asp residues comprising a negatively charged patch on the surface of the wtSH3. This undesirable behavior led to incorrect predictions regarding NMR measurablesparticularly, backbone and sidechain 15N relaxation rates in the flexible arginine-rich tail. After much experimentation with different MD setups and force fields, we choose to address this problem by slightly modifying the Lennard-Jones parameters of the pairwise N η−O ε/δ interactions associated with Arg-to-Glu/Asp salt bridges. In this connection, we noted a number of prior studies where Lennard-Jones potentials were selectively adjusted to improve the agreement between MD simulations and experiments (see paper 112). Very recently, this tactic has been used to improve the modeling of solvation for disordered proteins and peptides.19,35−37 In our study, a modest increase of 3% in the LJ equilibrium distance for Nη−Oε/δ pairs led to a successful MD model of the Sos:wtSH3 complex in good agreement with the experimental data. It is clear, however, that any such empirical correction requires careful cross-validation. We suggest that the Sos·dmSH3 complex examined in this paper offers a good opportunity to validate our approach. Indeed, as will be shown below, Sos·dmSH3 is a fuzzy complex, which is far more fluid and more dependent on electrostatic interactions than Sos:wtSH3. It therefore provides an independent test for the previously introduced force field correction. In pursuing the validation agenda, we first recorded a 5.4 μs trajectory of Sos·dmSH3 using the conventional (uncorrected) Amber ff99SB*-ILDN force field.13−15 The results from this MD trajectory, termed MD(xray,λ = 1.00), are presented in Figure 4. In this figure, panel A illustrates the positions of the peptide on the surface of the protein as sampled during the simulation. Gold spheres represent the center of mass of the hydrophobic N-terminal portion of Sos (PPPVPP), whereas blue spheres represent the center of mass of the positively charged Cterminal segment (RRRR). Clearly, the arginine-rich tail of Sos remains localized in the vicinity of the negatively charged patch on the surface of the dmSH3 (blue spheres cluster around RT and n-Src loops which are painted red in the graph). However, it samples a wider region of space than in the complex with wtSH3 (see paper 112). The N-terminus, on the other hand, is strongly delocalized. The typical binding pose is the one where the peptide is anchored to the dmSH3 via the C-terminus, while N-terminus sticks out into the solvent. This is a direct consequence of the Y186L/W169F mutations, which indeed disrupt the binding of the N-terminus as intended. Figure 4. Position of the peptide Sos on the surface of dmSH3 according to the data from the 5.4 μs simulation MD(xray, λ = 1.00). The gold and blue spheres represent the center of mass of the peptide N- and C-termini, respectively. The data are sampled with a 2 ns step, so that the number of spheres is approximately the same as in the corresponding graph of paper 1.12 The arrows and numbers in panel A indicate the two gold-sphere clusters which correspond to the distinctive binding poses discussed in the text. Six charged residues at the dmSH3 binding interface are colored red (D147, E148, E149, D150, E166, E167) and six adjacent hydrophobic residues are colored green (F141, F143, F169, P183, P185, and L186). Shown are the side view (A) and the top view (B) of the complex. Figure 5. Peptide−protein hydrogen bonds (A) and salt bridges (B) in the MD(xray, λ = 1.00) trajectory. In generating panel A, all hydrogen bonds that originate at the backbone of Sos and are populated at the level of at least 3% have been included. In generating panel B, we have taken into consideration all possible combinations between the 4 Arg residues from the Sos peptide and 6 Glu/Asp residues from the complementary acidic patch on the surface of the dmSH3 domain. Upon closer inspection we observe two distinct binding poses with respect to positioning of the N-terminus. One of them resembles the native-like conformation found in the Sos:wtSH3 (the corresponding cluster is marked with number 1 in the plot). This state originates from the initial X-ray-based coordinates (i.e., the coordinates of Sos:wtSH3 where the two point mutations have been introduced in silico). It survives during the first 0.4 μs of the trajectory, then disappears and never makes an appearance again. We note that the native-like binding pose remains a part of the dynamic equilibrium in the Sos·dmSH3 model but represents a local minimum which no longer dominates the thermodynamic landscape (cf. the illustration in Figure 1B). The second distinctive binding pose (cluster 2 in the plot) features residue P5 of the peptide packed against the aromatic ring of F169, but from the “wrong” side, i.e., the peptide is shifted to the left compared to the wild-type complex. In this pose, the segment P1−V4 is mostly projected into solvent, although it also makes on and off contacts with the hydrophobic surface of dmSH3. Additional information can be found in Figure 5, where the colored traces represent time evolution of the peptide-toprotein hydrogen bonds and salt bridges in the MD(xray, λ = 1.00) trajectory. Initially, the pattern of salt bridges resembles the one that was previously observed in the corresponding Sos:wtSH3 simulation (see paper 112). At around 0.4 μs this pattern changes; then at around 1.9 μs it changes again. At 2.7 μs the initial pattern is restored and persists, with some variations, to the end of the trajectory. Note, however, that during this time interval the system bears little resemblance to the wild-type complex; instead, the peptide assumes pose 2 which has been described above. The emergence of peptide-toprotein hydrogen bonds during this time interval, see Figure 5A, indicates that the conformational state of the peptide continues to evolve. Figure 5 also gives one an idea about the stability of salt bridges in the conventional (unaltered) MD simulation of Sos· dmSH3. The most prominent salt bridges, such as the nativelike bridge R7-D150, have a lifetime on the order of 1 μs. We have used the data in Figure 5 to calculate temporal correlation functions representing salt bridge formation/dissolution. These results (not shown) confirm the above rough estimate. To further characterize peptide binding, we have calculated the surface-to-surface distance between Sos and dmSH3 throughout the course of the MD(xray, λ = 1.00) trajectory (see Figure S3, panel A). The results indicate that the peptide remains bound to the protein over the entire duration of the simulation. The single fly off event occurs at 0.5 μs and lasts only 140 ps before the peptide rebinds at the same site. The fraction of the free Sos in the trajectory is, therefore, a miniscule 0.003%. This is far below what one may expect given the effective concentration of the protein and the peptide in this simulation, 8.1 mM, and the experimentally determined dissociation constant, ∼200 μM. To summarize, in the MD(xray, λ = 1.00) trajectory Sos remains strongly attached to dmSH3 at all times, primarily via the salt-bridge interactions. The pattern of the peptide-toprotein salt bridges is variable, with significant restructuring occurring on the time scale of microseconds. There are also adventitious nonpolar contacts that affect the binding poses of the peptide. One such contact, involving P5 and F169, is described above. We have also observed other such contacts in the additional λ = 1.00 simulations (not shown). The Nterminus of the peptide is typically projected into solvent and samples many different orientations. However, owing to its PPPVPP sequence, the N-terminal segment tends to form a Figure 6. Comparison of the experimental (red) and simulated MD(xray, λ = 1.00) (blue) 15N relaxation parameters and 1H,15N binding shifts for Sos peptide in the complex with dmSH3. Shown are the data from the backbone amide sites in the five nonproline residues in Sos as well as the data from four arginine side-chain 1Hε,15Nε sites. The 1Hε,15Nε signals are poorly resolved in the spectra and could not been assigned; therefore, the corresponding experimental data have been averaged. To facilitate the comparison with experiment, the respective simulated data have also been averaged. The sample conditions were 0.5 mM 15N-labeled Sos, 2 mM unlabeled dmSH3; all experimental values have been corrected to subtract out the contribution from 11% of the free peptide (see the Materials and Methods). The plotted data, therefore, pertain directly to the encounter complex Sos·dmSH3. PPII helix and shows only limited internal flexibility compared to a nonproline peptide. Finally, we have used the MD(xray, λ = 1.00) data to simulate 15N relaxation parameters as well as 1HN and 15N binding shifts of Sos. As already pointed out, the peptide remains bound to the protein during the simulation; therefore, the standard protocol to calculate 15N relaxation rates is appropriate. Briefly, the overall tumbling of the Sos·dmSH3 complex is first factored out and then reintroduced using the experimentally determined τrot = 7.0 ns. Other details of the computational protocol are described in the Materials and Methods. The inspection of Figure 6 makes it immediately obvious that MD data do not agree well with the experiment. In particular, this is evident for 15N transverse relaxation rates R2 and ηxy (panels B, F, and D), which are sensitive to slower motional modes (e.g., the overall molecular tumbling in the case of a rigid complex). The cross-correlated relaxation rate ηxy is especially relevant, since this parameter is free of exchange broadening and thus offers a clean measure of nanosecond dynamics. In light of these results, we can explain the origin of the discrepancy observed in Figure 6. The simulated data (blue circles in Figure 6, panels B, F, and D) suggest that Sos attaches itself to dmSH3 mainly through salt-bridge interactions involving its C-terminal arginines. In particular, residues R8 and R9 become essentially rigidly attached to the surface of the protein. The lifetime of the corresponding salt bridges, as well as backbone-to-sidechain hydrogen bonds, is longer than the tumbling time of Sos· dmSH3 (see Figure 5). As a result, R8 and R9 relax as if they were a part of the rigid peptide−protein complex, with the simulated backbone 15N R2 rates reaching 8 s−1. This is much higher than the experimentally measured R2 rates, which do not exceed 5 s−1 (red circles in Figure 6B). We suggest that the problem is caused by Arg-to-Asp/Glu salt bridges that proved to be excessively strong in the conventional MD(xray, λ = 1.00) simulation.12 The overstabilized salt bridges restrict and slow down the motion of the Sos peptide, which in turn leads to overestimated R2 and ηxy rates. This explanation is consistent with the 1HN,15N NOE data, which are sensitive to fast (picosecond) motional modes. The results suggest that fast motions are also overrestrained in the MD simulation, see Figure 6C,G. The same holds true for the side-chain 15Nε R1 rates, Figure 6E. As for the backbone data, one should bear in mind that dependence of R1 on motional correlation time is Λ-shaped; as a consequence, the simulated and experimental values happen to be similar in this case (see Figure 6A). The above notion is also supported by the binding shift data (see Figure 6, panels H and I). As already pointed out, the MD(xray, λ = 1.00) simulation features several relatively stable binding poses. This contributes to large predicted binding shifts Δδ in residues R7, R8, and R9. In reality, the system experiences a high degree of dynamic averaging which greatly reduces Δδ shifts. The emerging picture is similar to what has been previously found in the study of Sos:wtSH3 complex. In both cases MD Figure 7. Position of the peptide Sos on the surface of dmSH3 according to the data from the 5.2 μs simulation MD(xray, λ = 1.03). The plotting conventions are the same as in Figure 4; arrows and numbers in panel A identify the gold-sphere clusters, as described in the text. Figure 8. Peptide−protein hydrogen bonds (A) and salt bridges (B) in the MD(xray, λ = 1.03) trajectory. Plotting conventions are the same as in Figure 5. The interval of time between 4.7 and 5.0 μs corresponds to the binding pose 3, as described in the text; hydrogen bonds characteristic of this pose are not shown in the graph since they do not meet the 3% threshold. simulations overestimate the strength of the Arg-to-Asp/Glu salt bridges, leading to a misrepresentation of peptide dynamics. Therefore, we resort here to the same strategy as previously developed for the wild-type complexnamely, use the empirical force-field correction to selectively reduce the strength of the arginine salt bridges. It is worth pointing out that the two systems investigated in our two-part study, Sos:wtSH3 and Sos·dmSH3, are significantly different. Suffice it to mention that in Sos:wtSH3 the issue with salt bridges mainly affects the dynamic behavior of the two terminal residues, R9 and R10, whereas in Sos·dmSH3 the dynamics of the entire peptide is affected. This gives us reason to think that the proposed force-field correction has a sufficiently general significance.12,18 Simulation of Sos·dmSH3 Using Amber ff99SB*-ILDN with λ = 1.03 Correction: Comparison with Experimental Data. To address the problem of excessively strong Arg-toAsp/Glu salt bridges, we have introduced the empirical force field correction whereupon the Lennard-Jones equilibrium distances for Nη, Oε, and Nη, Oδ atom pairs is increased by 3%. Note that this correction is highly selective as it only influences the said salt bridges but not any other interactions involving Figure 9. Comparison of the experimental (red) and simulated MD(xray, λ = 1.03) (blue) 15N relaxation parameters and 1H,15N binding shifts for Sos peptide in complex with dmSH3. Plotting conventions are the same as in Figure 6. these residues, such as interactions with water. The resulting (minimally) adapted version of the force field has been used to record 5.2 μs trajectory of the Sos·dmSH3 complex starting from the X-ray-based coordinates. Figure 7 illustrates the positioning of the Sos peptide over the surface of the dmSH3 domain in the MD(xray, λ = 1.03) simulation. The peptide is strongly delocalized in this simulation. Nevertheless, it remains clear that the binding is dominated by electrostatic interactions between the positively charged C-terminus of Sos and the negatively charged RT and n-Src loops in dmSH3 (cf. massive blue clusters in the upper left part of Figure 7A). The N-terminus of Sos for the most part is projected into solvent, as evident from wide scatter of the gold spheres in the graph. At the same time, several distinct binding poses, similar to those that have been discussed in the previous section, can be identified here as well. One cluster, which is indicated by an arrow and number 1 in Figure 7A, corresponds to the native-like peptide binding pose, as inherited from the initial set of coordinates. Similar to what has been described above, this pose survives until ∼0.4 μs when N-terminus loses connection with the peptide-binding hydrophobic patch. Later, this pose is revisited at 3.3 μs and again at 5.0 μs; both visits are brief, lasting for ∼0.1 μs. Another cluster, which is marked by number 2 in the plot, corresponds to a relatively stable pose that persists from 4.1 to 4.7 μs in the trajectory. In this pose the C-terminal segment of the peptide is attached to its usual binding spot, whereas the Nterminus is locked at the site far away from the nominal binding surface. Specifically, the peptide is held in place by the bifurcated hydrogen bond between the N-terminus of P1 and side-chain carboxylic sites of D174 and E176, as well as hydrogen bonds between the side-chain of R10 and backbone amide sites of N144 and D147. Transient salt bridges between R7, R9 on one hand and E148, E149 on the other hand are also a part of this binding pose (see Figure 8). Interestingly, the Cterminus of the peptide tends to form a short 310-helix in this binding pose. Finally, the cluster marked with number 3 corresponds to the time interval from 4.7 to 5.0 μs when the peptide becomes coordinated to the opposite face of dmSH3. During this time period, Sos is attached to the protein mainly via transient salt bridges and hydrogen bonds formed by R8, R9, and R10 with several C-terminal residues of the protein as well as hydrogen bonds between P5 or P6 and the side chain of Y190. These and other similar interactions are highly dynamic; none of them persist for too long, but together they help to keep Sos localized at this adventitious site for about 0.3 μs. A hydrophobic contact between two of Sos prolines and Y136 is also a part of this binding pose. Figure 8 sheds additional light on the behavior of the simulated encounter complex. Clearly, the lifetime of the intermolecular Arg-to-Asp/Glu salt bridges in MD(xray, λ = 1.03) trajectory is much shorter than what has been previously observed. The corresponding temporal correlation functions (not shown) demonstrate that the more prominent salt bridges remain stable for no more than tens of nanoseconds. This is an order-of-magnitude drop compared to the unmodified MD(xray, λ = 1.00) simulationsalt bridges are dissolved and re-formed much faster in the new trajectory. Intermolecular hydrogen bonds remain stable for more extended intervals of time, although without the continuing stabilizing effect of the salt bridges they are also less visible than before, cf. Figures 5 and 8. Finally, the analysis of the simulated 15N relaxation rates and binding shifts indicates that MD(xray, λ = 1.03) trajectory is in excellent agreement with the experimental data, Figure 9. This result is in marked contrast with the prior unaltered MD simulation. Briefly, in the previously discussed MD(xray, λ = 1.00) trajectory, Sos was too strongly bound to dmSH3 via multiple salt bridges. As a consequence, the simulated relaxation rates in Sos were not very different from what can be observed in a rigid protein−peptide complex with the tumbling time of 7 ns. Furthermore, the predicted binding shifts were also similar to what may be expected from a tight complex, i.e., on the order of 1 ppm. In contrast, in the new MD(xray, λ = 1.03) simulation the peptide is significantly less constrained by the salt bridges. While Sos remains attracted to the negatively charged patch on the surface of the dmSH3, it is at the same time much more dynamic, rapidly “dancing” on the surface of the protein. Its effective motional correlation time is considerably less than 7 ns, resulting in lower R2 and ηxy rates. The contacts between Sos and dmSH3 are highly variable and labile; as a consequence, binding shifts are effectively averaged out, with predicted Δδ values on the order of 0.1 ppm. The emerging picture is that of a fuzzy complex, where the dynamic status of Sos is somewhere in between the free and fully bound peptide. This is further illustrated in Figure 11 and visualized in the videos, which are provided as a part of Supporting Information. As discussed above, certain distinctive binding poses are observed in the MD(xray, λ = 1.03) trajectory; however, these poses are less prominent than in the previous MD(xray, λ = 1.00) simulation and also more dynamic and conformationally diverse. At this point it is appropriate to discuss the dissociation and rebinding of Sos in the MD(xray, λ = 1.03) simulation. As it turns out, Sos spends 2.9% of the entire simulation time separated from dmSH3 (see Figure S3, panel B). This is 1000 times more than in the MD(xray, λ = 1.00) trajectory, where the peptide binding proved to be excessively tight. The fact that Sos occasionally loses connection with dmSH3 is relevant for MD-based 15N relaxation rate calculations. Recall that MD trajectories have been processed such as to eliminate the overall rotation of the dmSH3. The overall rotation was subsequently reintroduced into dipolar correlation functions g(τ) through the multiplicative factor exp(−τ/τrot). This procedure is valid when the peptide maintains direct contact with the protein, i.e., as long as the Sos·dmSH3 complex tumbles as a single entity. However, this approach is strictly speaking inapplicable in relation to those time intervals when Sos becomes separated from dmSH3. We notice, though, that the peptide spends only a small fraction of time in separation from the protein, so that it should not have any appreciable impact on the calculated 15N relaxation rates. Furthermore, during the time intervals when Sos becomes detached from dmSH3, it moves as a free peptide, i.e., shows rapid internal dynamics on the time scale less than 1 ns. On the other hand, the described postprocessing scheme deals with a longer time scale associated with the overall tumbling time τrot. Therefore, the postprocessing treatment has little impact on g(τ) insofar as free Sos is concerned. In essence, our MD-based calculations of the 15N relaxation rates are soundthey simply ignore a small correction arising from 2.9% of the free peptide. To further test the validity of this approach, we have implemented an alternative procedure to calculate Sos 15N relaxation rates. Specifically, we have extracted dipolar correlation functions g(τ) directly from the MD trajectories, without any prior processing. As a next step, we redefined the time axis of g(τ), multiplying it by a factor 2.7. This step is intended to correct for the intrinsically low viscosity of the TIP3P water model.33,34 The results (not shown) were in reasonably good agreement with the experimental data, although not as good as those in Figure 9. More specifically, the alternative strategy works well for those relaxation rates that significantly depend on the overall rotational motion, i.e., R2 and ηxy. At the same time, small but appreciable discrepancies are observed for relaxation parameters that are sensitive to fast internal motions, i.e., R1 and NOE. In any event, it is important that MD(xray, λ = 1.03) data show significantly better agreement with the experimental results than MD(xray, λ = 1.00) also using this alternative processing strategy. Next, we turn to the discussion of Kd and kOFF as can be estimated on the basis of the MD trajectory. As has already been mentioned, the fraction of free Sos in the MD(xray, λ = 1.03) trajectory is 2.9%. The effective concentration of Sos and dmSH3 in the simulation is 8.1 mM. Hence the dissociation constant for the in silico Sos·dmSH3 complex is 7 μM, which is somewhat lower than the experimental value of ∼200 μM. The discrepancy is, however, relatively modest, corresponding to just 2 kcal/mol difference in binding free energy. This is comparable to the accuracy of the most sophisticated MDbased methods for calculation of protein−peptide interactions.38,39 The reality is probably more complicated than the above simple argument. Consider, for example, Sos peptide that lost contact with dmSH3 (in a sense of Figure S3), yet some of its Arg side chains remain within the electrostatic capture radius of the Glu/Asp patch. In this situation Sos is likely to rebind. Therefore, this arrangement should be viewed as a part of Sos· dmSH3 ensemble, rather than a combination of a free Sos and dmSH3. Generally, we invoke two different approaches to Kd determination: (i) the experiment-based approach using chemical shift changes at the individual amide sites and (ii) the MD-based approach using simple distance threshold. We do not necessarily expect these two approaches to produce identical results. Currently, we can only state that MD(xray, λ = 1.03) simulation is qualitatively consistent with the experimental data with respect to both Kd and kOFF estimations. In addition to the data presented in Figure 9H-J, it is also instructive to consider the reverse titration data, i.e., changes in amide chemical shifts of dmSH3 originating from its interaction with Sos. Because Sos does not contain any aromatic rings and does not form any stable hydrogen bonds with dmSH3, protein resonances exhibit only modest binding shifts, see Figure S2. These shifts are too small to allow for any meaningful comparison with MD/SHIFTX240 predictions. However, they can be compared with the corresponding wtSH3 data, see paper 1.12 As it turns out, in the case of Sos:wtSH3 the biggest shifts are observed on residues V184, W169, Q168, N146, and D147. The average Δδ value for these six protein residues is 169 Hz (with proton and nitrogen shifts combined according to Δδ = (ΔδH2 + ΔδN2)1/2). All these residues are compactly located around the side chain of W169. On the other hand, in the case of Sos·dmSH3, the residues showing biggest binding shifts are V184, V187, E149, and Q168; the average Δδ value for those residues is 63 Hz. These sites are more broadly distributed on the surface of the protein, consistent with the picture of the fuzzy complex (see Figure S2). Finally, we will consider an additional set of experimental parameters, namely, temperature coefficients associated with 1 N H chemical shifts. The corresponding peptide data for Sos:wtSH3, Sos·dmSH3, and free Sos are summarized in Table 2. In general, 1HN temperature coefficients (TCs) more positive than −4.5 ppb/K indicate that the amide is involved in a hydrogen bond, whereas TCs more negative than −4.5 ppb/K correspond to solvent-exposed amides.27 Table 2. Experimental 1HN Temperature Coefficients for Residues in Sos Peptide (ppb/K) Sos:wtSH3 Sos·dmSH3 free Sos V4 R7 R8 R9 R10 −7.45 −9.32 −9.56 −5.87 −6.64 −6.75 −3.88 −5.46 −7.13 −5.15 −6.14 −6.12 −6.54 −6.47 −6.49 First, it is instructive to compare the data from Sos:wtSH3 with those from free Sos peptide (top and bottom rows, respectively). The TC value of R10 proves to be the same for the free and the bound peptide, −6.5 ppb/K. Indeed, in the wild-type complex this terminal residue belongs to the dangling peptide tailit is fully solvated, has little interaction with the protein, and behaves similar to free peptide (see paper 112). All other residues, however, show significantly different TC values. In particular, for R8 the TC increases from −7.13 to −3.88 ppb/K upon formation of the complex. This increase is symptomatic of the intermolecular hydrogen bond formed by the amide group of R8 and the side-chain carboxyl of E166. A substantial increase is also observed for R7, R9, and especially V4, reflecting significant contacts that these residues make with wtSH3. Next, let us compare the data from Sos·dmSH3 with those from free Sos peptide (middle and bottom rows, respectively). As it turns out, the TC values in those two samples are virtually indistinguishable, with the exception of residue R8, which shows an appreciable difference. This result is consistent with the previous titration data indicating that the binding shifts in Sos·dmSH3 are very small except for R8 (see Figure 3). It is also consistent with the MD-based picture of fuzzy complex wherein Sos is highly dynamic and forms only transient contacts with dmSH3. The MD simulation, however, lacks predictive power to identify the origin of the effect observed in R8. It can be suggested that, on average, backbone amide site of R8 tends to interact more strongly with dmSH3. Alternatively, one may suggest that proximity to dmSH3 affects the conformational distribution of Sos (e.g., promotes intramolecular hydrogen bonding) and thus brings about the observed effect in R8. This possibility is discussed in more detail in the next section. Simulation of Sos Binding to dmSH3: Convergence of the MD Results. Figure 10 shows the initial scene from MD(rand, λ = 1.03) simulation. Displayed in this plot is the simulated Sos−dmSH3 pair (central molecules) along with the periodic images of dmSH3. Only those images that lie in the xz plane are shown; the grid spacing along the y axis is the same. The space between protein molecules is filled with TIP3P solvent. One can appreciate that there is a significant room left for the peptide to fully explore a range of conformations. Any periodic boundary artifacts, e.g., a situation where Sos interacts with the “master copy” of dmSH3 and at the same time with one of its periodic images, should be minimal. The initial coordinates shown in Figure 10 have been generated as described in paper 1.12 Briefly, a random conformation of Sos is extracted from 2.1 μs trajectory of the free peptide and then placed at the distance of ∼15 Å above the surface of dmSH3. For the trajectory at hand, MD(rand, λ = 1.03), the starting placement of the peptide is such that Nterminal residues P1 and P2 initially obstruct the access to the RT loop of dmSH3. After 1 ns, these residues move sideways and the first salt bridge, R8-D147, is formed. Subsequently, the Figure 10. Initial placement of Sos and dmSH3 in the 6.0 μs MD(rand, λ = 1.03) trajectory recorded in the explicit TIP3P solvent. Shown is the 3 × 3 square arrangement, containing the central (simulated) protein molecule along with its multiple periodic images. The spacing between the protein molecules in the y dimension (not shown) is the same as in x and z dimensions, with the minimal surface-to-surface separation of 40.6 Å. The MD(xray, λ = 1.00) and MD(xray, λ = 1.03) models have similar dimensions (cf. Table 1). Sos peptide behaves very similar to the MD(xray, λ = 1.03) simulation: it spends most of the time “dancing” over the surface of dmSH3, sometimes finding relatively stable binding poses and occasionally experiencing brief fly off episodes. The 15N relaxation rates and binding shifts calculated from MD(rand, λ = 1.03) simulation are shown in Figure S4. The results are very similar to those shown in Figure 9, suggesting that our simulations of Sos·dmSH3 achieve reasonably good convergence. The most significant discrepancy in the MD(rand, λ = 1.03) results involves 15N binding shift of residue R7. This shift is predicted to be 1.41 ppm, which is much larger than the experimentally determined value of −0.06 ppm, see Figure S4. On the other hand, the previously discussed MD(xray, λ = 1.03) trajectory predicts the shift of 0.03 ppm, which is in good agreement with the experiment. What is the cause of this inconsistency? Considering the conformational ensemble of free Sos (2.1 μs trajectory), we find that 41% of the time the peptide features a distinctive secondary-structure element. Specifically, Sos forms a short helical turn that starts at residue P5 and ends at R9. In some of the MD frames this turn can be classified as a 310 helix. As such, it is defined by P5-R8 and, to a lesser extent, P6-R9 backbone hydrogen bond; it is further stabilized by the backbone-to-sidechain hydrogen bond involving V4 and R9. In other MD frames, this turn resembles the so-called Némethy helix,41 featuring the bifurcated hydrogen bond P5-R8/R9. This helical element is relatively stable, existing on the time scale of hundreds of nanoseconds. If we now turn to the MD(xray, λ = 1.03) trajectory we will find the same helical turn, which happens to be populated at the same level, 41%. This turn persists for extended periods of time: for instance, it is present from 4.1 to 4.7 μs when the peptide assumes binding pose 2 (see discussion above). Because the content of this helical turn is the same in the free-Sos trajectory and in MD(xray, λ = 1.03), it is not reflected in the calculated binding shiftsthe corresponding contribution is subtracted out. On the other hand, the content of the helical turn in the MD(rand, λ = 1.03) trajectory proves to be a mere 11%. When calculating binding shifts, ΔδSos = δbound − δfree Sos Sos , the effect of helicity is no longer canceled. It is well-known that 15N experiences a significant upfield shift when the peptide is in helical conformation.42 As a consequence, the calculation produces larger than expected value of Δδ. The effect is clearly of statistical naturethe two supposedly equivalent trajectories of Sos·dmSH3 exhibit a different proportion of the helical turn, 41% vs 11%. Of note, an additional random-start trajectory recorded in our study shows 34% fraction of the helical turn. Apparently, this particular element of the transient secondary structure is relatively longlived; the current ∼5 μs trajectories are insufficient to generate adequate statistics in this regard. Not surprisingly, chemical shifts prove to be especially sensitive to the lack of convergence concerning residual secondary structure. A similar argument can be made with respect to binding poses. As discussed above, the lifetime of certain poses approaches 1 μs. Clearly, longer trajectories are needed to obtain a proper statistical sampling of these poses. Nevertheless, the current set of simulations proves to be successful in capturing the remarkable dynamic nature of Sos· dmSH3, as summarily illustrated in the final section. ■ CONCLUSION A combination of experimental NMR with MD simulations has proven to be a highly successful strategy to study dynamics of folded proteins. In broad terms, NMR data provide validation for MD simulations; in turn, MD simulations offer valuable information that cannot be easily gleaned from the NMR data. This strategy holds even greater promise for intrinsically disordered proteins. However, it also faces bigger challenges. First, extra-large water boxes are needed to accommodate the fully extended conformations of disordered proteins and extralong trajectories are needed to sample the vast conformational space available to IDPs. Until recently this was beyond the reach of computer hardware. Second, currently available force fields have been designed and optimized with folded proteins in mind and have never been systematically tested for IDPs. Therefore, it remains unclear whether such simulations can be trusted. In this situation, much effort has been invested in developing alternative approaches beyond Molecular Dynamics, e.g., ensemble models of IDPs.43−45 This situation has started to change in recent years. A number of MD simulations have been conducted for disordered peptides and proteins.2,8,9,46−50 Multiple recent studies undertook a systematic comparison of different force fields in the context of protein disorder and folding.37,51−60 In many cases a good agreement with experiment has been noted, although certain systematic deviations have also been found. In the last 2 years, this led to a number of force-field developments targeted specifically at IDPs.19,61−63 Given that experimental information on IDPs is intrinsically limited (in particular, highresolution structural data are unavailable), it can be anticipated that in future all-atom MD simulations will greatly enhance the potential of NMR experiments in this area. Along the same lines, in this two-part study of Sos:wtSH3 and Sos·dmSH3 complexes, which both feature elements of conformational disorder, we have identified certain shortcomings in the MD force fields. As it turns out, the strength of the Arg-to-Asp/Glu salt bridges is not correctly reproduced (typically overestimated) in the most widely used conventional force fields. It is not accidental that this finding has been made in the context of conformational disorder. Indeed, the problem can be easily overlooked in the simulations of folded proteins. It is known that the contributions of salt bridges to protein stability are context-dependent and generally small.64−66 Therefore, it is difficult to use thermodynamic data to unmask the issue with force-field parametrization of salt bridges.67 The problem tends to escape detection, unless specifically targeted.68−70 In contrast, in our study of the fuzzy encounter complex we have found that the strength of the salt bridges directly determines the dynamic behavior of the peptide. Consequently the problem with force field parametrization became readily apparent. We predict that λ = 1.03 correction should also work well for other protein systems because of its highly selective character it targets Arg-to-Asp/Glu salt bridges without perturbing any other interactions in the system. It is expected that this result should stimulate further efforts in the area of force field development. Along these lines, D. Cerutti is currently collaborating with K. Debiec and others to develop a pairspecific Lennard−Jones matrix which should properly balance the strength of polar interactions.71 In a long term, efforts in this direction should greatly benefit from the advent of polarizable force fields.72 Figure 11. Summary of 15N backbone relaxation data from Sos as obtained in the study of the stable Sos:wtSH3 complex (paper 112), dynamic Sos· dmSH3 complex and free Sos (this paper). Shown on the left are the snapshots from the three respective MD trajectories that symbolically represent the relevant conformational ensembles. Shown on the right is the comparison between the experimental and simulated (λ = 1.03) relaxation rates from all three systems. Shown in Figure 11 is the summary of 15N backbone relaxation data from three distinct species: Sos tightly bound to the wtSH3, Sos loosely associated with dmSH3, and free Sos. In terms of dynamic behavior, Sos:wtSH3 and free Sos represent the two extremes. The former is bulky and tumbles slowly, resulting in high R2 and ηxy rates (green profiles in the plot). The latter is small and flexible, resulting in low R2 and ηxy rates (black dashed-line profiles). As it turns out, the engineered complex of Sos·dmSH3 falls between these two extremes. In this fuzzy complex the peptide engages in a sort of a “break dance” on the surface of the protein. The characteristic time scale of peptide motion in this case is shorter than in the conventional Sos:wtSH3 complex but longer than in free Sos. This remarkable dynamic situation produces the medium-range R2 and ηxy rates (magenta profiles). After we have implemented a targeted force-field correction to remedy the problem with Arg-to-Asp/Glu salt bridges, our MD simulations nicely reproduce the same hierarchy of dynamics species as observed experimentally (see Figure 11). This is particularly encouraging given that our data set contains two distinctly different systems. Sos:wtSH3 is an essentially rigid complex, where only the arginine-rich peptide tail is dynamic and sensitive to the strength of salt-bridge interactions. On the other hand, Sos·dmSH3 is an essentially fluid complex, where dynamic status of the entire peptide critically depends on the strength of the salt-bridge interactions. In both cases we have achieved good agreement between the MD simulations and experiment. We believe that we have thereby produced a realistic MD model of the electrostatic encounter complex, which is guided and supported by the experimental NMR data. Importantly, this model relies on backbone and side-chain 15N relaxation data. Relaxation data are indispensable to elucidate the dynamic status of disordered systems. Of note, relaxation data can be rigorously modeled only on the basis of bona f ide MD simulations but not on the basis of static conformational ensembles.8,45 In this sense our study is in line with the emerging trend in the field which progresses from ensemble models of protein disorder to full-fledged MD simulations. ■ ASSOCIATED CONTENT S Supporting Information * The Supporting Information is available free of charge on the ACS Publications website at DOI: 10.1021/acs.biochem.5b01283. Table with rotational diffusion parameters for wild-type and double-mutant c-Crk N-SH3 in the free form and in the presence of Sos peptide; superposition of the 1 N 15 H - N HSQC spectra of the wtSH3 and dmSH3 domain; 1HN-15N HSQC titration for dmSH3 domain binding Sos and the corresponding chemical shift perturbation map; figure illustrating the separation between Sos and dmSH3 in the MD(xray, λ = 1.00) and MD(xray, λ = 1.03) trajectories; figure presenting the comparison of the experimental and simulated 15N relaxation parameters and 1H,15N binding shifts in Sos· dmSH3 using the MD(rand, λ = 1.03) trajectory (PDF) Animations of the MD(xray, λ = 1.00) simulation of Sos· dmSH3: the entire trajectory with the duration 5.4 μs (sampling step 10 ns) and a fragment of the trajectory with the duration 0.1 μs (sampling step 100 ps); animations of the MD(xray, λ = 1.03) simulation of D Sos·dmSH3: the entire trajectory with the duration 5.2 μs (sampling step 10 ns) and a fragment of the trajectory with the duration 0.1 μs (sampling step 100 ps) (MPG) ■ REFERENCES (1) Fawzi, N. L., Phillips, A. H., Ruscio, J. Z., Doucleff, M., Wemmer, D. E., and Head-Gordon, T. (2008) Structure and dynamics of the Aβ21−30 peptide from the interplay of NMR experiments and molecular simulations. J. Am. Chem. Soc. 130, 6145−6158. (2) Sgourakis, N. G., Merced-Serrano, M., Boutsidis, C., Drineas, P., Du, Z. M., Wang, C. Y., and Garcia, A. E. (2011) Atomic-level characterization of the ensemble of the Aβ1−42 monomer in water using unbiased Molecular Dynamics simulations and spectral algorithms. J. Mol. Biol. 405, 570−583. (3) Rembert, K. B., Paterova, J., Heyda, J., Hilty, C., Jungwirth, P., and Cremer, P. S. (2012) Molecular mechanisms of ion-specific effects on proteins. J. Am. Chem. Soc. 134, 10039−10046. (4) Knott, M., and Best, R. B. (2012) A preformed binding interface in the unbound ensemble of an intrinsically disordered protein: evidence from molecular simulations. PLoS Comput. Biol. 8, e1002605. (5) Allison, J. R., Rivers, R. C., Christodoulou, J. C., Vendruscolo, M., and Dobson, C. M. (2014) A relationship between the transient structure in the monomeric state and the aggregation propensities of α-synuclein and β-synuclein. Biochemistry 53, 7170−7183. (6) Karp, J. M., Erylimaz, E., and Cowburn, D. (2015) Correlation of chemical shifts predicted by molecular dynamics simulations for partially disordered proteins. J. Biomol. NMR 61, 35−45. (7) Ball, K. A., Phillips, A. H., Nerenberg, P. S., Fawzi, N. L., Wemmer, D. E., and Head-Gordon, T. (2011) Homogeneous and heterogeneous tertiary structure ensembles of amyloid-β peptides. Biochemistry 50, 7612−7628. (8) Xue, Y., and Skrynnikov, N. R. (2011) Motion of a disordered polypeptide chain as studied by paramagnetic relaxation enhancements, 15N relaxation, and Molecular Dynamics simulations: how fast is segmental diffusion in denatured ubiquitin? J. Am. Chem. Soc. 133, 14614−14628. (9) Lindorff-Larsen, K., Trbovic, N., Maragakis, P., Piana, S., and Shaw, D. E. (2012) Structure and dynamics of an unfolded protein examined by Molecular Dynamics simulation. J. Am. Chem. Soc. 134, 3787−3791. (10) Robustelli, P., Trbovic, N., Friesner, R. A., and Palmer, A. G. (2013) Conformational dynamics of the partially disordered yeast transcription factor GCN4. J. Chem. Theory Comput. 9, 5190−5200. (11) Mantsyzov, A. B., Maltsev, A. S., Ying, J. F., Shen, Y., Hummer, G., and Bax, A. (2014) A maximum entropy approach to the study of residue-specific backbone angle distributions in alpha-synuclein, an intrinsically disordered protein. Protein Sci. 23, 1275−1290. (12) Xue, Y., Yuwen, T. R., Zhu, F. Q., and Skrynnikov, N. R. (2014) Role of electrostatic interactions in binding of peptides and intrinsically disordered proteins to their folded targets. 1. NMR and MD characterization of the complex between the c-Crk N-SH3 domain and the peptide Sos. Biochemistry 53, 6473−6495. (13) Hornak, V., Abel, R., Okur, A., Strockbine, B., Roitberg, A., and Simmerling, C. (2006) Comparison of multiple Amber force fields and development of improved protein backbone parameters. Proteins: Struct., Funct., Genet. 65, 712−725. (14) Best, R. B., and Hummer, G. (2009) Optimized Molecular Dynamics force fields applied to the helix-coil transition of polypeptides. J. Phys. Chem. B 113, 9004−9015. (15) Lindorff-Larsen, K., Piana, S., Palmo, K., Maragakis, P., Klepeis, J. L., Dror, R. O., and Shaw, D. E. (2010) Improved side-chain torsion potentials for the Amber ff99SB protein force field. Proteins: Struct., Funct., Genet. 78, 1950−1958. (16) Wu, X. D., Knudsen, B., Feller, S. M., Zheng, J., Sali, A., Cowburn, D., Hanafusa, H., and Kuriyan, J. (1995) Structural basis for the specific interaction of lysine-containing proline-rich peptides with the N-terminal SH3 domain of c-Crk. Structure 3, 215−226. (17) Andrews, C. T., and Elcock, A. H. (2013) Molecular Dynamics simulations of highly crowded amino acid solutions: comparisons of eight different force field combinations with experiment and with each other. J. Chem. Theory Comput. 9, 4585−4602. (18) Debiec, K. T., Gronenborn, A. M., and Chong, L. T. (2014) Evaluating the strength of salt bridges: a comparison of current biomolecular force fields. J. Phys. Chem. B 118, 6561−6569. (19) Best, R. B., Zheng, W. W., and Mittal, J. (2014) Balanced protein-water interactions improve properties of disordered proteins and non-specific protein association. J. Chem. Theory Comput. 10, 5113−5124. (20) Baker, C. M., and Best, R. B. (2014) Insights into the binding of intrinsically disordered proteins from molecular dynamics simulation. WIREs Comput. Mol. Sci. 4, 182−198. (21) Sekhar, A., and Kay, L. E. (2013) NMR paves the way for atomic level descriptions of sparsely populated, transiently formed biomolecular conformers. Proc. Natl. Acad. Sci. U. S. A. 110, 12867− 12874. (22) Demers, J. P., and Mittermaier, A. (2009) Binding mechanism of an SH3 domain studied by NMR and ITC. J. Am. Chem. Soc. 131, 4355−4367. (23) Lovell, S. C., Word, J. M., Richardson, J. S., and Richardson, D. C. (2000) The penultimate rotamer library. Proteins: Struct., Funct., Genet. 40, 389−408. (24) Guerois, R., Nielsen, J. E., and Serrano, L. (2002) Predicting changes in the stability of proteins and protein complexes: A study of more than 1000 mutations. J. Mol. Biol. 320, 369−387. (25) Schymkowitz, J., Borg, J., Stricher, F., Nys, R., Rousseau, F., and Serrano, L. (2005) The FoldX web server: an online force field. Nucleic Acids Res. 33, W382−W388. (26) Muhandiram, D. R., and Kay, L. E. (1994) Gradient-enhanced triple-resonance three-dimensional NMR experiments with improved sensitivity. J. Magn. Reson., Ser. B 103, 203−216. (27) Baxter, N. J., and Williamson, M. P. (1997) Temperature dependence of 1H chemical shifts in proteins. J. Biomol. NMR 9, 359− 369. (28) Kroenke, C. D., Loria, J. P., Lee, L. K., Rance, M., and Palmer, A. G. (1998) Longitudinal and transverse 1H-15N dipolar 15N chemical shift anisotropy relaxation interference: Unambiguous determination of rotational diffusion tensors and chemical exchange effects in biological macromolecules. J. Am. Chem. Soc. 120, 7905−7915. (29) Jorgensen, W. L., Chandrasekhar, J., Madura, J. D., Impey, R. W., and Klein, M. L. (1983) Comparison of simple potential functions for simulating liquid water. J. Chem. Phys. 79, 926−935. (30) Case, D. A., Darden, T. A., Cheatham, T. E., Simmerling, C. L., Wang, J., Duke, R. E., Luo, R., Walker, R. C., Zhang, W., Merz, K. M., Roberts, B., Wang, B., Hayik, S., Roitberg, A., Seabra, G., Kolossváry, I., Wong, K. F., Paesani, F., Vanicek, J., Liu, J., Wu, X., Brozell, S. R., Steinbrecher, T., Gohlke, H., Cai, Q., Ye, X., Wang, J., Hsieh, M.-J., Cui, G., Roe, D. R., Mathews, D. H., Seetin, M. G., Sagui, C., Babin, V., Luchko, T., Gusarov, S., Kovalenko, A., Kollman, P. A. (2010) AMBER 11, University of California, San Francisco, CA. (31) Case, D. A., Babin, V., Berryman, J. T., Betz, R. M., Cai, Q., Cerutti, D. S., Cheatham, T. E., Darden, T. A., Duke, R. E., Gohlke, H., Goetz, A. V., Gusarov, S., Homeyer, N., Janowski, P., Kaus, J., Kolossváry, I., Kovalenko, A., Lee, T. S., LeGrand, S., Luchko, T., Luo, R., Madej, B., Merz, K. M., Paesani, F., Roe, D. R., Roitberg, A., Sagui, C., Salomon-Ferrer, R., Seabra, G., Simmerling, C. L., Smith, W., Swails, J., Walker, R. C., Wang, J., Wolf, R. M., Wu, X., Kollman, P. A. (2014) AMBER 14, University of California, San Francisco, CA. (32) Humphrey, W., Dalke, A., and Schulten, K. (1996) VMD: visual molecular dynamics. J. Mol. Graphics 14, 33−38. (33) van der Spoel, D., van Maaren, P. J., and Berendsen, H. J. C. (1998) A systematic study of water models for molecular simulation: derivation of water models optimized for use with a reaction field. J. Chem. Phys. 108, 10220−10230. (34) Gonzalez, M. A., and Abascal, J. L. F. (2010) The shear viscosity of rigid water models. J. Chem. Phys. 132, 096101. (35) Nerenberg, P. S., and Head-Gordon, T. (2011) Optimizing protein-solvent force fields to reproduce intrinsic conformational preferences of model peptides. J. Chem. Theory Comput. 7, 1220−1230. (36) Nerenberg, P. S., Jo, B., So, C., Tripathy, A., and Head-Gordon, T. (2012) Optimizing solute-water van der Waals interactions to reproduce solvation free energies. J. Phys. Chem. B 116, 4524−4534. (37) Henriques, J., Cragnell, C., and Skepo, M. (2015) Molecular Dynamics simulations of intrinsically disordered proteins: force field evaluation and comparison with experiment. J. Chem. Theory Comput. 11, 3420−3431. (38) Gumbart, J. C., Roux, B., and Chipot, C. (2013) Standard binding free energies from computer simulations: what is the best strategy? J. Chem. Theory Comput. 9, 794−802. (39) Gumbart, J. C., Roux, B., and Chipot, C. (2013) Efficient determination of protein-protein standard binding free energies from first principles. J. Chem. Theory Comput. 9, 3789−3798. (40) Han, B., Liu, Y. F., Ginzinger, S. W., and Wishart, D. S. (2011) SHIFTX2: significantly improved protein chemical shift prediction. J. Biomol. NMR 50, 43−57. (41) Kuster, D. J., Liu, C., Fang, Z., Ponder, J. W., and Marshall, G. R. (2015) High-resolution crystal structures of protein helices reconciled with three-centered hydrogen bonds and multipole electrostatics. PLoS One 10, e0123146. (42) Wishart, D. S., Sykes, B. D., and Richards, F. M. (1991) Relationship between Nuclear Magnetic Resonance chemical shift and protein secondary structure. J. Mol. Biol. 222, 311−333. (43) Marsh, J. A., and Forman-Kay, J. D. (2012) Ensemble modeling of protein disordered states: Experimental restraint contributions and validation. Proteins: Struct., Funct., Genet. 80, 556−572. (44) Schneider, R., Huang, J.-R., Yao, M., Communie, G., Ozenne, V., Mollica, L., Salmon, L., Jensen, M. R., and Blackledge, M. (2012) Towards a robust description of intrinsic protein disorder using nuclear magnetic resonance spectroscopy. Mol. BioSys. 8, 58−68. (45) Ball, K. A., Wemmer, D. E., and Head-Gordon, T. (2014) Comparison of structure determination methods for intrinsically disordered amyloid-β peptides. J. Phys. Chem. B 118, 6405−6416. (46) Mao, A. H., Crick, S. L., Vitalis, A., Chicoine, C. L., and Pappu, R. V. (2010) Net charge per residue modulates conformational ensembles of intrinsically disordered proteins. Proc. Natl. Acad. Sci. U. S. A. 107, 8183−8188. (47) Sethi, A., Tian, J. H., Vu, D. M., and Gnanakaran, S. (2012) Identification of minimally interacting modules in an intrinsically disordered protein. Biophys. J. 103, 748−757. (48) Stanley, N., Esteban-Martin, S., and De Fabritiis, G. (2014) Kinetic modulation of a disordered protein domain by phosphorylation. Nat. Commun. 5, 5272. (49) Levine, Z. A., Larini, L., LaPointe, N. E., Feinstein, S. C., and Shea, J.-E. (2015) Regulation and aggregation of intrinsically disordered peptides. Proc. Natl. Acad. Sci. U. S. A. 112, 2758−2763. (50) Yedvabny, E., Nerenberg, P. S., So, C., and Head-Gordon, T. (2015) Disordered structural ensembles of vasopressin and oxytocin and their mutants. J. Phys. Chem. B 119, 896−905. (51) Best, R. B., and Mittal, J. (2011) Free-energy landscape of the GB1 hairpin in all-atom explicit solvent simulations with different force fields: Similarities and differences. Proteins: Struct., Funct., Genet. 79, 1318−1328. (52) Cino, E. A., Karttunen, M., and Choy, W. Y. (2012) Effects of molecular crowding on the dynamics of intrinsically disordered proteins. PLoS One 7, e49876. (53) Vymetal, J., and Vondrasek, J. (2013) Critical assessment of current force fields. Short peptide test case. J. Chem. Theory Comput. 9, 441−451. (54) Gerben, S. R., Lemkul, J. A., Brown, A. M., and Bevan, D. R. (2014) Comparing atomistic molecular mechanics force fields for a difficult target: a case study on the Alzheimer’s amyloid β-peptide. J. Biomol. Struct. Dyn. 32, 1817−1832. (55) Piana, S., Klepeis, J. L., and Shaw, D. E. (2014) Assessing the accuracy of physical models used in protein-folding simulations: quantitative evidence from long molecular dynamics simulations. Curr. Opin. Struct. Biol. 24, 98−105. (56) Skinner, J. J., Yu, W., Gichana, E. K., Baxa, M. C., Hinshaw, J. R., Freed, K. F., and Sosnick, T. R. (2014) Benchmarking all-atom simulations using hydrogen exchange. Proc. Natl. Acad. Sci. U. S. A. 111, 15975−15980. (57) Fluitt, A. M., and de Pablo, J. J. (2015) An analysis of biomolecular force fields for simulations of polyglutamine in solution. Biophys. J. 109, 1009−1018. (58) Hoffmann, K. Q., McGovern, M., Chiu, C.-C., and de Pablo, J. J. (2015) Secondary structure of rat and human amylin across force fields. PLoS One 10, e0134091. (59) Palazzesi, F., Prakash, M. K., Bonomi, M., and Barducci, A. (2015) Accuracy of current all-atom force-fields in modeling protein disordered states. J. Chem. Theory Comput. 11, 2−7. (60) Rauscher, S., Gapsys, V., Gajda, M. J., Zweckstetter, M., De Groot, B. L., and Grubmüller, H. (2015) Structural ensembles of intrinsically disordered proteins depend strongly on force field: a comparison to experiment. J. Chem. Theory Comput. 11, 5513−5524. (61) Wang, W., Ye, W., Jiang, C., Luo, R., and Chen, H. F. (2014) New force field on modeling intrinsically disordered proteins. Chem. Biol. Drug Des. 84, 253−269. (62) Piana, S., Donchev, A. G., Robustelli, P., and Shaw, D. E. (2015) Water dispersion interactions strongly influence simulated structural properties of disordered protein states. J. Phys. Chem. B 119, 5113− 5123. (63) Mercadante, D., Milles, S., Fuertes, G., Svergun, D. I., Lemke, E. A., and Grater, F. (2015) Kirkwood-Buff approach rescues overcollapse of a disordered protein in canonical protein force fields. J. Phys. Chem. B 119, 7975−7984. (64) Kumar, S., and Nussinov, R. (2002) Close-range electrostatic interactions in proteins. ChemBioChem 3, 604−617. (65) Makhatadze, G. I., Loladze, V. V., Ermolenko, D. N., Chen, X. F., and Thomas, S. T. (2003) Contribution of surface salt bridges to protein stability: Guidelines for protein engineering. J. Mol. Biol. 327, 1135−1148. (66) Lee, C. W., Wang, H. J., Hwang, J. K., and Tseng, C. P. (2014) Protein thermal stability enhancement by designing salt bridges: a combined computational and experimental study. PLoS One 9, e112751. (67) Geney, R., Layten, M., Gomperts, R., Hornak, V., and Simmerling, C. (2006) Investigation of salt bridge stability in a Generalized Born solvent model. J. Chem. Theory Comput. 2, 115−127. (68) Kaas, Q., Aumelas, A., Kubo, S., Chino, N., Kobayashi, Y., and Chiche, L. (2002) The [Lys−2-Arg−1-des(17−21)]-endothelin-1 peptide retains the specific Arg−1-Asp8 salt bridge but reveals discrepancies between NMR data and molecular dynamics simulations. Biochemistry 41, 11099−11108. (69) Gruia, A. D., Fischer, S., and Smith, J. C. (2003) Molecular dynamics simulation reveals a surface salt bridge forming a kinetic trap in unfolding of truncated Staphylococcal nuclease. Proteins: Struct., Funct., Genet. 50, 507−515. (70) Leontyev, I. V., and Stuchebrukhov, A. A. (2010) Electronic continuum model for molecular dynamics simulations of biological molecules. J. Chem. Theory Comput. 6, 1498−1508. (71) Cerutti, D. S., Swope, W. C., Rice, J. E., and Case, D. A. (2014) ff14ipq: a self-consistent force field for condensed-phase simulations of proteins. J. Chem. Theory Comput. 10, 4515−4534. (72) Huang, J., Lopes, P. E. M., Roux, B., and MacKerell, A. D. (2014) Recent advances in polarizable force fields for macromolecules: microsecond simulations of proteins using the classical Drude oscillator model. J. Phys. Chem. Lett. 5, 3144−3150.