Article pubs.acs.org/JPCB
How Hydrophilic Proteins Form Nonspecific Complexes Ozlem Ulucan and Volkhard Helms* Center for Bioinformatics, Saarland University, 66041 Saarbrücken, Germany S Supporting Information *
ABSTRACT: In the crowded environment of cells, proteins frequently encounter other proteins in many possible orientations. Most of these encounters are short-lived because the physicochemical properties of the two binding patches do not match. However, even for protein pairs that bind tightly, it is not an easy task to find the correct binding site on the partner protein and align with it. So far not well understood is the source of interaction specificity that favors a small set of specific “native” interactions over the multitude of alternative orientations. We used molecular dynamics simulations to study nonspecific complexes formed by barnase and barstar, cytochrome c and cytochrome c peroxidase, and the complex of the Nterminal domain of enzyme I with the histidine-containing phosphocarrier. We found that spontaneously forming nonspecific encounters involve interaction interfaces smaller than those of the specific complexes and are attracted by shorter-range direct interactions between the proteins.
■
INTRODUCTION Nonspecific encounters of diffusing proteins frequently occur in living cells but have received far less attention than specific protein complexes. Interestingly, nonspecific encounters may even play important roles in the formation of specifically bound conformations of protein complexes. For example, site-directed mutagenesis1 and Brownian dynamics (BD) simulations2 of the barnase−barstar (BN−BS) pair showed that perturbations of the charge distribution outside the native binding surface can modulate the association rate. Utilizing BD simulations, Northrup et al.3 revealed that favorable electrostatic interactions facilitate long-lived nonspecific encounters that, at a later stage, convert to the reactive specific complex. On the basis of paramagnetic nuclear magnetic resonance spectroscopy, Volkov et al.4 estimated that complexes formed by the proteins cytochrome c (CC) and cytochrome c peroxidase (CCP) spend ∼30% of their lifetimes in dynamic encounter states. Using the paramagnetic relaxation enhancement (PRE) technique, Tang et al.5 demonstrated that even a brief, imperfect collision can mediate the formation of the specific complex between the Nterminal domain of enzyme I (EIN) and the histidinecontaining phosphocarrier (HPr). Their data suggest that protein surfaces not involved directly in the specific binding interface may facilitate the assembly of the functional complex. On the basis of data from PRE and replica exchange simulations, Hummer and co-workers6 reported that distinct nonspecific complexes exist for this system that account for ∼10% of the relative population. They also pointed out that besides accelerating the binding kinetics, nonspecific complexes also play a role in protein function as alternative binding modes. Despite these advances, our understanding of the structural and energetic aspects of nonspecific complex formation is quite limited. © XXXX American Chemical Society
We have recently employed restrained molecular dynamics (MD) simulations with umbrella potentials to characterize the one-dimensional association free energy pathways between specific complexes of the hydrophilic protein pairs of BN and BS, CC and CCP, and EIN and HPr and their unbound states.7 For comparison, we present here the findings from analogous MD simulations of nonspecific complexes formed by BN and BS, CC and CCP, and EIN and HPr.
■
METHODS The coordinates for the three specific protein−protein complexes were retrieved from the Protein Data Bank (PDB):8 barnase−barstar (PDB entry 1BRS9), cytochrome c−cytochrome c peroxidase (PDB entry 2PCC10), and the amino-terminal domain of enzyme I and the histidinecontaining phosphocarrier protein (PDB entry 3EZB11). The first set of simulations consisted of plain unbiased MD simulations that started from different distances that aimed to monitor the spontaneous assembly of nonspecific complexes (for details, see below). Umbrella sampling simulations were then conducted to compute the potential of mean force for the dissociation process starting from the identified encounter complexes. Trajectories were analyzed as described previously.7 All simulations were conducted using the Gromacs package version 4.5.412 and all-atom Amber force field FF99SB-ILDN.13 The force field parameters of the heme groups in the cytochrome c−cytochrome c peroxidase system were taken from our previous study.7 All MD simulation parameters (time step, use of LINCS, etc.), the potential of mean force (PMF) Received: June 18, 2015 Revised: July 28, 2015
A
DOI: 10.1021/acs.jpcb.5b05831 J. Phys. Chem. B XXXX, XXX, XXX−XXX
Article
The Journal of Physical Chemistry B Table 1. Selected Global and Interface Parameters of Nonspecific Complexes Observed during MD Simulationsa complex BN−BS
CC−CCP
EIN−HPr
native NSC−r1−x4 NSC−r1−y2 NSC−r1−y3 NSC−r1−y4 native NSC−r1−y4 NSC−r2−y2 native NSC−r1−y3 NSC−r1−y4 NSC−r2−x3 NSC−r2−y2 NSC−r2−y3
lifetimeb (ns) 20 92+ 31 50 27 51+ 92+ 73+ 50 77+ 45+
area of binding interface (Å2)
no. of interface residues
± ± ± ± ± ± ± ± ± ± ± ± ± ±
16:14 10:11 15:14 11:15 11:11 13:10 10:10 11:10 33:24 15:13 13:12 13:12 17:16 10:12
655 313 429 410 372 449 336 340 639 425 414 486 633 441
30 51 31 43 46 42 36 33 58 128 64 53 57 51
no. of H-bonds across interface 11.6 2.7 2.4 4.6 4.8 6.9 4.2 7.2 4.7 6.6 4.6 5.6 6.5 9.3
± ± ± ± ± ± ± ± ± ± ± ± ± ±
1.8 1.8 1.4 1.5 1.5 1.9 1.1 1.4 1.6 2.9 1.8 1.5 1.6 2.5
no. of salt bridges across interface
total charge of interface I:II (e)
± ± ± ± ± ± ± ± ± ± ± ± ± ±
+3:−4 0:−2 0:−3 0:−5 +2:−3 +5:−2 −1:+1 −2:+4 −5:+4 −3:+2 +1:+1 −1:+3 −2:+3 −3:0
4.0 1.0 0.6 0.5 0.9 3.6 1.8 1.3 3.8 2.2 0.3 2.2 3.4 0.5
0.6 0.9 0.4 0.7 0.7 1.0 0.8 0.7 1.2 1.1 0.5 0.8 0.7 0.5
a The number of interface residues calculated from the trajectories is based on cutoffs for the distance and frequency of occurrence (see Methods). The interface area and the number of H-bonds across the interfaces were calculated using the Gromacs utilities g_sas and g_hbond. The number of salt bridges was calculated using VMD. The nonspecific complexes that were selected as starting structures for the PMF simulations are shown in bold. bThe plus sign indicates that the proteins remained bound until the end of the 100 ns simulation after forming the initial contact.
simulation trajectories using visual inspection to identify those trajectories where at least a single contact is formed between the two proteins. Afterward, we assessed the binding events observed in these trajectories regardless of whether these events resulted in nonspecific complexes based on the following three criteria: (1) The contact should last at least 20 ns. The precise contact time is termed the lifetime of the nonspecific complex. (2) During this lifetime, each protein should contribute at least 10 residues to the interface. Also, each of these residues should participate at the interface in at least 80% of the MD snapshots sampled during the lifetime of the complex. This number was chosen for the following reason: in control simulations for the respective specific complexes, the interface residues of the X-ray complex maintained the contacts with the interface residues of the partner protein in at least 80% of the MD snapshots sampled during a 100 ns simulation. We defined those residues as interface residues that contained at least one single heavy atom within 0.5 nm of the partner protein. (3) The resulting complex should not be nativelike (for this, we superimposed the native and putative nonspecific complexes using the larger protein partner as a reference; then all putative nonspecific complexes that showed overlapping residues with the native complex were discarded). For each system, we selected those two nonspecific complexes for further analysis that bear either the longest lifetime or the largest area of the contact interface compared to the other nonspecific complexes formed. For each of them, the coordinates from the first snapshot of the ≥20 ns lifetime were taken as initial coordinates for umbrella sampling simulations.
protocol, and the methods of postanalysis to retrieve thermodynamic data from the umbrella potential simulations were equivalent to those described in our previous study.7 The first set of simulations consisted of plain MD simulations starting with a pair of unbound proteins to obtain the nonspecific complexes. The initial configurations for this set of simulations were generated by displacing one of the proteins to an interfacial distance of 1.0−2.5 nm and rotating it by 90°, 180°, and 270° about the x or y axis (see Table S1). If the rotation resulted in a clash of the two proteins, the interfacial distance was increased 0.5 or 1.5 nm. For each system, seven independent molecular dynamics simulations were conducted, including the unrotated configuration. Each simulation was repeated once after assignment of different random initial velocities. In this way, 14 simulations were conducted for each system. We also conducted 100 ns MD simulations of the specific bound complexes to characterize the behavior of interfacial residues over time. Each initial configuration was placed in a cubic box of TIP3P14 water of suitable dimensions so that at least 1.4 nm of solvent was around the two proteins in each direction. Thus, the size of the simulation boxes varied according to the initial separation of the two proteins. To mimic physiological conditions, a NaCl concentration of 100 mM was realized, including neutralizing counterions. Following an initial energy minimization of 1000 steps of steepest descent, each system was equilibrated in two steps whereby the heavy atoms of the proteins were restrained. The first step involved 100 ps of MD in the NVT ensemble, maintaining the temperature at 310 K. Protein and nonprotein atoms were coupled separately to temperature baths using Berendsen’s weak coupling algorithm.15 Subsequently, 100 ps NPT equilibration was performed, keeping the pressure at 1 bar also using Berendsen’s weak coupling method.15 During data collection, the Nose− Hoover thermostat16,17 was combined with the Parrrinello− Rahman barostat18 to regulate temperature and pressure, respectively. Data collection was completed by a 100 ns conventional MD simulation in the absence of any restraints. The definition of nonspecific complexes used here is arbitrary and applied for practical purposes. We analyzed the 100 ns MD
■
RESULTS AND DISCUSSION Interfacial Characteristics of Nonspecific Complexes. Of 14 MD simulations 100 ns in length for each protein pair, we observed binding events in 11, 7, and 13 MD simulations of the BN−BS, CC−CCP, and EIN−HPr systems, respectively. Of the bound conformations formed during these simulations, we consider four BN−BS complexes, two CC−CCP complexes, and five EIN−HPr complexes as “different enough” from the respective specific complex that these conformations likely represent nonspecific complexes. B
DOI: 10.1021/acs.jpcb.5b05831 J. Phys. Chem. B XXXX, XXX, XXX−XXX
Article
The Journal of Physical Chemistry B
Figure 1. Cartoon representations (left) of the nonspecific complexes selected for further analysis together with their specific counterparts. The larger component of each pair (BN, CCP, and EIN are colored blue) was utilized to superimpose the specific and two nonspecific complexes. The smaller component (BS, CC, and HPr) of the specific complex is colored red, and its orientations in the nonspecific complexes are colored orange and silver. Potentials of mean force (right) for the complexes shown in the left panel. The plain vertical lines represent the cutoff that separates bound and unbound regions. The vertical dashed lines mark the positions of the bound states.
number of H-bonds considerably larger than the number in nonspecific complexes. For some other systems, we observed nonspecific complexes that have a larger mean number of Hbonds (the NSC−r2−y2 nonspecific complex of the CC−CCP system and all nonspecific complexes of the EIN−HPr system except NSC−r1−y4) than their native counterparts. All nonspecific complexes contain fewer salt bridges across their interface than their specific counterparts. On the other hand, all interfaces of nonspecific complexes bear at least three charged residues, although some of them carry zero overall charges. This indicates the importance of electrostatic interactions, which were previously reported to be important in the nonspecific binding process.5,6,20 One-Dimensional Free Energy Profile and Standard Free Energy of Binding. Figure 1 depicts the onedimensional free energy profiles for the associations of the nonspecific complexes computed from umbrella potential simulations. The PMF curves exhibit overall a similar behavior. The physical separations where the PMF curves start to flatten differ among the systems and are generally shorter for nonspecific complexes than for the respective specific complexes. The picture is somewhat different for the EIN− HPr system. For this system, the PMF curve flattens out at physical separations of 1.3 nm (NSC−r1−y3) and 1.5 nm (NSC−r2−y2), which is similar to the case for the specific
Table 1 lists general and interface characteristics of the nonspecific complexes formed. The interface properties reported for the native complexes (the area of the binding interface and the number of H-bonds and salt bridges) in this work slightly differ from the values we reported previously,7 for the following two reasons: (1) The values reported in this work were calculated over a trajectory of 100 ns rather than a single frame, and (2) different programs were used to perform the calculations. H-Bonds were determined on the basis of cutoffs for the donor−hydrogen−acceptor angle (
