ARTICLE pubs.acs.org/JPCB
Exploring the Thermodynamic Landscape, Kinetics, and Structural Evolution of a Protein Conformational Transition with a Microscopic Double-Well Model Zai-Zhi Lai,† Qiang Lu,† and Jin Wang*,†,‡,§,^ †
Department of Chemistry, ‡Department of Physics, and §Department of Applied Mathematics & Statistics, State University of New York at Stony Brook, Stony Brook, New York 11794, United States ^ State Key Laboratory of Electroanalytical Chemistry, Changchun Institute of Applied Chemistry, Chinese Academy of Sciences, Changchun, Jilin 130021, People’s Republic of China ABSTRACT: Functional conformational transition in the glutamine-binding protein (GlnBP) is known to be the key to bind and transfer ligand glutamine. Here, we developed a structure-based double-well model to investigate the thermodynamic and kinetic natures of the GlnBP conformational transition. We uncovered the underlying free-energy landscape of the conformational transition with different temperatures. The analysis shows that below the melting temperature, two basins of attractions emerge, corresponding to the open state and the closed state of the protein. We explored the kinetic property of the conformational switch through the mean and distribution of the first passage time as well as the autocorrelation function. The kinetics implies the complexity and the hierarchical structure of the underlying energy landscape. We built the contact maps of the structures to probe the structural evolution of the conformational transition. Finally, the φ values of the residues were calculated to identify the important residues (hot spots) of the transition state.
’ INTRODUCTION Transitions between different conformational states are the key for the biological function of many biomolecules. Protein flexibility and plasticity allow the protein to bind ligands, form oligomers, and perform mechanical work. Current experiments are able to characterize the structures of biomolecule using X-ray diffraction and NMR spectroscopy. The dynamical properties of biomolecules can be probed by spectroscopic techniques.1,2 However, obtaining the detailed pictures of biomolecular dynamics, especially for the global structural conformational transition information, by experimental techniques remains a longstanding endeavor. Computational simulation provides another approach to study the structure, motion, and function of biomolecules. It potentially provides full time-dependent structural information on biomolecules. Yet, most processes of interest occur on time scales (microsecond to second) inaccessible to standard all-atom molecular dynamics simulation. One approach to overcome the problem of long time scales in simulation is to use a simplified model. We propose developing a structure-based model to achieve a long-time molecular dynamics simulation by describing the protein interactions in a coarse-grained (CG) way at the residue level in which the water molecules are not explicitly included and native interactions are preferred. Energy landscape theory36 provides a theoretical framework for understanding many biochemical processes, such as protein conformational transitions and protein folding. The topology of the energy landscape affects the dynamic behavior of biochemical processes and influences functional and stability properties. r 2011 American Chemical Society
The relationship between the dynamic, stability, and functional behavior and the structure of the protein can be quantitatively understood by the detailed description of the free-energy landscape. One important concept of the standard energy landscape theory is the funnel-like shape of the landscape, which controls the kinetic and thermal properties of the folding process. For the process of the conformational transition, the multiple native structures of the protein require the multibasin shape of the underlying energy landscape. A Go-like model,7 which was proposed to emphasize the importance of the native structure, can be applied to describe to topology of the energy landscape. This model considers that the global shape of the free-energy landscape seems mostly determined by the native structure and applies the Go-type interactions, that is, the attraction interactions are assigned to the pair of residues that interact in the native structure, and repulsive interactions are endowed to the other contacts. The model was expected to provide useful information about the topology of the energy landscape. In this article, a new interaction was considered in the standard structure-based model to simulate the large conformational changes and to study the dynamic properties controlled by the underlying multiple-basin energy landscape of proteins. The new interaction forms a twowell potential, which includes the two minima between two contact residues, one for the open and one for the closed structure, and also a barrier separating these two minima. Received: November 12, 2010 Revised: January 31, 2011 Published: March 22, 2011 4147
dx.doi.org/10.1021/jp110845u | J. Phys. Chem. B 2011, 115, 4147–4159
The Journal of Physical Chemistry B
ARTICLE
In this article, we used glutamine-binding protein (GlnBP) as the model protein to simulate the conformational transitions. Two reference structures were supplied by X-ray crystallography,8,9 as shown in Figure 1a and b. For the ligand-free open structure of GlnBP, the PDB code is 1GGG; for the ligand-bound closed structure of GlnBP, the PDB code is 1WDN. GlnBP is a small periplasmic binding protein in Escherichia coli. It contains a single polypeptide chain of 226 residues. The tertiary structure of GlnBP consists of 35% Rhelix and 37% β-sheet. GlnBP is composed of two similar globular domains. The large domain includes two separate peptide segments, that is, residues 184 and residues 186 226, and the small domain includes residues 90180. These two domains are linked by two peptide hinges, which are residues 8589 and residues 181185.
’ MODEL The potential energy U for a given protein conformation Γ in the simulation is N1
UðΓ, Γ1 , Γ2 Þ ¼
∑
Kb ðbi bonds N2 þ
Kθ ðθi θ1i Þ2 ðθi θ2i Þ2
φ þ φ2i φ φ2i 2 cos φi 1i cos 1i 2 2 !12 nonnat native C þ E þ Unat ðrij Þ ð1Þ rij ji jj > 3 ji jj > 3
N3
þ
∑
angles
b0i Þ2
∑ Kφ dihedrals
∑
∑
where Γ1(Γ2) represents the native open(closed) state. The first three terms describe the deviation of bond lengths, bond angles, and torsion angles away from equilibrium values. The nonbonding interaction potential can be implemented by the sum of the nonnative interaction and the native interaction, which are given by the fourth term and fifth terms. For the bonded term, due to the change of the distance between the open state and closed state being very small, a single harmonic shape potential was applied to calculate the bonded interaction. For the angles term and dihedrals term, the change between the open state and closed state is significant; therefore, we extended the standard forms to the two-well forms to calculate the corresponding interactions. The explicit representation of the native interactions is given below10
8 k > r1 > > E ZðrÞðZðrÞ aÞ with ZðrÞ ¼ if r < r1 > 1 > r > > > n > > Y ðrÞ > > Y ðrÞ ¼ ðr rh Þ2 ðrh r1 Þ2n > > n 2 > > þ E CY ðrÞ with 4nðE1 þ E2 Þ if r1 e r < rh 2 > > 2n C¼ > > > r2 rh > > > 2ðm 1Þ < B ¼ E1 mðr2 rh Þ > Eh ðm 1Þðr2 rh Þ2 > > Y ðrÞ h1 h1 ¼ > > B with if rh e r < r2 mðEh þ E2 Þ > > Y ðrÞm þ h2 > > 2m > E ðm 1Þðr r Þ > 2 2 h > > h2 ¼ > > E h þ E2 > " > # > 12 10 > > r2 r2 > > > if r2 e r E2 5 6 > : r r
(
ð2Þ
where the potential energy is a function of the interacting distance r, between two residues; r1 represents the position of the minimum of the first potential well, r2 is the position of the minimum of the second potential well, and rh is the location of the maximum barrier of the potential energy. To obtain the smooth and continuous curve for the shape of the potential energy, the parameters m = 5, k = 8, and n = 1 are used in this paper. The schematic representation of the potential is shown in Figure 2. The first well corresponds to the closed structure, and the second well corresponds to the open structure. It is worth mentioning that this two-well model is microscopic. It represents the interaction between each two residues. This is different from the previous studies,1113 and we can term it as a macroscopic two-well model. The interactions between two residues have a single well, while the resulting macroscopic states can have two free-energy wells with preferred open and closed structures. Here, our model can be called a microscopic two-well model. The microscopic interactions between two residues are two wells, resulting macroscopically in two or more wells including the open and closed states. We consider the distance between and those in the closed state two residues in the open state ropen ij . When ropen = rclose , a single-well potential is applied to the rclose ij ij ij rclose |< residues i and j, as shown in Figure 2a; when 0 < |ropen ij ij 2 Å, a two-well potential with a shallow barrier emerges to the residues i and j, and the form of this kind of two-well potential is very similar to the single well, as shown in Figure 2b; in the case of rclose | > 2 Å, the two-well potential emerges with a |ropen ij ij significant barrier, as shown in Figure 2c. The native contact pairs for the open and the closed states are derived from the CSU software.14 There are three parameters ɛ1, ɛ2, and ɛh that we can choose to determine the form of the two-well potential, as shown in Figure 2c. The ɛ1 represents the depth of the first well. This depth mainly controls the melting temperature of the protein. The ɛ2 is the depth of the second well. The difference between ɛ1 and ɛ2, that is, Δɛ = |ɛ1 ɛ2|, controls the relative stability of the two states. The last parameter ɛh is the energy barrier of the two states. In our model, the melting temperature is close to 355 K. Here, ɛ1 was set to 0.58 kcal/mol. The Δɛ and ɛh can be determined from experimental data. For convenience, we adjusted these two parameters so that the reasonable conformational switch happened between the open state and the closed state. After fixing the temperature at 310 K, we kept changing the value of Δɛ until the protein spent almost equal time in each state. Then, we gradually increased the value of ɛh to control the emergence of the transitions within the reasonable computation time. Here, we set Δɛ = 0.12 kcal/mol and ɛh = 0.35 kcal/mol.
’ METHODS AND RESULTS After fixing the parameters, we chose four temperatures, 275, 310, 336, and 355 K, to perform the simulations. For each of these four temperatures, 50 trajectories with 50 million steps for each trajectory were simulated using modified AMBER software in order to obtain reliable statistical results. Additionally, to observe the behaviors of the protein at different temperatures, we also performed 25 trajectories for temperatures of 280, 290, 300, 320, and 350 K. The protein began to melt at 355 K in our model. Figure 3 shows the typical trajectories at temperatures 275, 310, 336, and 355 K. Below the melting temperature, we can see the 4148
dx.doi.org/10.1021/jp110845u |J. Phys. Chem. B 2011, 115, 4147–4159
The Journal of Physical Chemistry B
ARTICLE
Figure 1. (a) The unbound-open structure of GlnBP. The region with the green color is the large domain of GlnBP, and the region with the yellow color is the small domain. The hinge region is blue. (b) The bound-closed structure of GlnBP. The green and yellow regions represent the large domain and the small domain of GlnBP, respectively. The hinge region is blue, and the ligand glutamine is marked in red.
reversible transition between two states. The conformational switch occurs very rapidly, without any intermediate state. In other words, the breaking of the contacts in the initial structure and formation of the contacts in the finial structure occur simultaneously. This observation is consistent with the previous study.12 The root mean square deviations (rmsds) of the two states exhibit that the protein in the open state has greater conformational fluctuation than that in the closed state. This implies that, with the ligand binding, the proteinligand complex may be more stable than the ligandunbound structure. Interestingly, from Figure 3ac, one can observe that the protein prefers to stay in the closed state at low temperature and biases to the open state at high temperature. This phenomenon can also be observed clearly in Figure 4a, which describes quantitatively the percent chance of staying in the open state or the closed state of the protein at different temperatures. At 275 K, the protein has insufficient energy, so that the kinetic transition between open and closed states occurs very slowly. The protein remains in the closed state for 80% of the simulation time. At 310 K, the protein spends almost the same time in both states. When the temperature reaches to 336 K, the conformational change occurs frequently, and the open state is dominant. The essential nature of the protein determines its thermodynamic behavior. In our microscopic two-well potential, the closed state has lower energy than the open state. At low temperature, the protein does not have enough energy to break the ligand-induced contacts. In other words, it is difficult for the
Figure 2. The shapes of the potentials. (a) Single-well potential. ɛ1 is the depth of the well. (b) Two-well potential with shallow barrier. ɛ1 is the depth of the first well; ɛ2 is the depth of the second well; ɛh is the height of the barrier. (c) Two-well potential with a significant barrier. ɛ1 is the depth of the first well; ɛ2 is the depth of the second well; ɛh is the height of the barrier. In this paper, |ɛ1| > |ɛ2|. The height of the barrier in (b) is half of that in (c).
protein to jump out of the first well ɛ1. Therefore, most of the time, the protein is mostly trapped in the first well at low temperature. As the temperature increases, the energy of the protein increases, and it has more chance to go over the barrier; therefore, the percentage of staying in the closed state decreases. At high temperature, the thermal energy of the protein can break all of the ligand-induced contacts quickly when they are formed. Additionally, from Figure 4b, we can know that at low temperature, the average resident time of the closed state is longer than that of the open state, and the 4149
dx.doi.org/10.1021/jp110845u |J. Phys. Chem. B 2011, 115, 4147–4159
The Journal of Physical Chemistry B
ARTICLE
Figure 3. Typical trajectories of the simulations at temperatures (a) 250, (b) 310, (c) 336, and (d) 355 K. The red color represents the rmsd of the closed state, and the green color represents the rmsd of the open state.
average residence times of both states decrease when the temperature increases. The two-dimensional free-energy profiles, as shown in Figure 5, were constructed by the trajectories at temperatures 275, 310, 336, and 355K. The free energy was considered as a function of rmsd1 and rmsd2 and obtained by F = log(P), where P is the statistical population obtained from all 50 trajectories for each temperature. In Figure 5ac, there are two local minima in each profile corresponding to the open state (right basin) and the closed state (left basin) of the protein. In Figure 5a, there are few conformations connecting the two basins, implying that the transition is slow at this low temperature. Also, the basin of the closed state is deeper than that of the open state, which suggests that the protein prefers to stay in the closed state at low temperature. In Figure 5b, the depths of the two minima are almost the same, and the height of the single free-energy barrier between the two minima equals 4.4kBT. One feature of this free-energy landscape is that the right basin corresponding to the open state is broader than the left basin corresponding to the closed state. Another interesting topological property is the steepness of the basins. In Figure 5b, the left basin is steeper than the right basin. This topological characteristic illustrates the following transition dynamics. When the conformational transition starts from the open state, it will undergo a relatively more extensive pathway to reach the transition state, which will then go downhill to the closed state rapidly on the energy surface. The distinct feature of the freeenergy profile in Figure 5d is that it has a long “valley”, except for the two basins, suggesting the melting state at this temperature. In contrast, the closed-state basin in Figure 5d is more shallow than the open-state basin. The conformational transition is important for biological function of GlnBP. The kinetics of the transition between the open state and closed state can be characterized by the distribution of the first passage time (FPT) of open states, P(To), where To is the dwell time of the open state, and closed states, P(Tc), where Tc is the closed residence time, respectively. In particular,
derived from the FPT distribution, the mean first passage times R ¥ of the open state and closed state, that is, ÆTo(c)æ 0 TP(To(c)) dT, are significant quantifiers of the transition kinetics. Specifically, they determine the mean opening(closing) rates ro(c) ÆTc(o)æ1. Furthermore, one may connect the distribution of the FPT to the other kinetic properties, such as the barrier of the transition, controlled by the underlying energy landscape. Figure 6 shows the first passage time distributions for the open state and closed state at temperature 275, 310, and 336 K. We fitted the distributions with Γ distributions, P(To(c)) TRo(c) exp(λTo(c)), which take the asymptotic form P(To(c)) exp(λo(c)To(c)) for large To(c), where the λo(c) are the exponents corresponding to the open state and closed state, respectively, and R is a real number. As the temperature increases, both λo and λc increase; λo is 4.5 104(275 K), 1.2 103(310 K), and 1.5 103(336 K), and λc equals 1.0 104(275 K), 9.5 104(310 K), and 6.7 103(336 K), as shown in Figure 6. The increase of the exponent λ suggests that the rate of the conformational switch increases. In addition, at the high temperature, the difference between the value of the mean FPT and the maximal point of the FPT distribution is much smaller than that at the low temperature (data not shown), which means that for the low temperature, the FPT distribution provides more of a detailed description of the kinetic properties of the conformational transition. The dynamical motion of the conformation switch can also be analyzed by calculating the autocorrelation function c(k), where k is the time interval separating the same variable (in rmsd relative to the native state (open state or closed state)) from one time to another. The formula that we used to calculate the c(k) is the following n1
cðkÞ ¼
∑ ðxt ÆxæÞðxt þ k ÆxæÞ
t¼1
nk
∑
t¼1
4150
ðxt ÆxæÞ
ð3Þ
2
dx.doi.org/10.1021/jp110845u |J. Phys. Chem. B 2011, 115, 4147–4159
The Journal of Physical Chemistry B
ARTICLE
Figure 4. (a) Percentages of staying in two states at different temperatures. The red line corresponds to the open state, and the green line corresponds to the closed state. (b) Average residence time of two states at different temperatures. The red line corresponds to the open state, and the green line corresponds to the closed state.
where xt is the value of the rmsd to the native state at time t and Æxæ is the average rmsd value of the whole trajectory. Using the above formula to calculate the c(k) needs a large n, the number of data. Therefore, we simulated long trajectories with 500 millions steps at 310 K (unless specified otherwise, the following results in this paper are all at a temperature of 310 K), which provides the statistical reliability to calculate the autocorrelation function c(k). Figure 7 shows the autocorrelation function in linear and semilogarithmic coordinates. In the semilogarithmic coordinates, the y value is logarithmic, and the x value remains linear. We see from the Figure 7a and c that the autocorrelation function decays with time interval k. We also can clearly observe from Figure 7b and d that the decay shows exponential behavior at the short time limit, k < 1500. In this region, the correlation function can be fitted quit well with a single-exponential function. The time scale
is 450 for the closed state, and is 451 for the open state, as shown in Figure 7a and c. When k > 1500, the correlation function cannot be fitted with a single-exponential function. The multiexponential function can fit the autocorrelation function well in the long time regime (represented in semilogarithmic coordinates), indicating that the kinetics of the whole conformation transition process is multiexponential. At 310 K, the average open-dwell time scale approaches 1200, and the average closedwell time scale approaches 1300, as shown in Figure 4b. From Figure 7d, we can see that for the open state, there are three correlation time scales, that is, 499, 1000, and 1206. Two of them are close to the average open-dwell time, and the other one is shorter than that, yet it is very close to the time scale of the single-exponential function in Figure 7c, implying that the intrabasin transition is faster than the interbasin transitions from 4151
dx.doi.org/10.1021/jp110845u |J. Phys. Chem. B 2011, 115, 4147–4159
The Journal of Physical Chemistry B
ARTICLE
Figure 6. Distributions of the first passage time of the different temperatures. Every x axis in (af) represents the number of the molecular simulation steps of staying in the open and closed states.
Figure 5. Free-energy profile of the conformational transition of GlnBP at temperatures (a) 250, (b) 310, (c) 336, and (d) 355 K. The x and y axes are the rmsds to the closed and open structures, respectively.
the open state to closed state. For the closed state, as shown in Figure 7b, the three correlation time scales are 347, 1304, and 1521. Two of them approach the average kinetic time, and the other one is close to the single-exponential function for the closed state, suggesting that the interbasin transition is the slow process that dominates the average kinetics. Because the kinetics can be considered as a good probe of the underlying energy
landscape, the multiexponential kinetics implies the complexity of the energy landscape. When a complex system with these characteristics of kinetics diffuses on the corresponding underlying energy landscape from one state to the other, it will depend on a series of subprocesses. A conformational switch can be recognized as a consequence of these many subprocesses, all of which make the transitions possible.15,16 To probe the structural evolution in the process of the conformational transition, the residue contact pairs of 11 different structures on the pathway from the open state to the closed state were analyzed by the residue contact map.17 First, we chose 11 points (states) on the pathway of the conformational transition on the free-energy landscape. The distance between any two adjacent points (states) is moderate, so that these points can cover the whole pathway. These 11 points included the open state, the closed state, and the transition state. These three states correspond to the 1st, 6th, and 11th structures of all 11 structures that we chose. Next, we calculated the probabilities of the contact pairs in different states. The probability means how closely a specific pair accesses to its native distances (open state or closed state). The formula used is Pkij = nkij/nk, where Pkij represents the probability of the pair ij in state k, nk is the number of conformations around the state k, and nkij is the weight of the pair ij at state k. The weight was calculated by the distances of the pair ij in the conformations around state k. Figures 811 show the contact maps and their corresponding structures. In these figures, the right column shows the contact 4152
dx.doi.org/10.1021/jp110845u |J. Phys. Chem. B 2011, 115, 4147–4159
The Journal of Physical Chemistry B
ARTICLE
Figure 7. Autocorrelation function distribution of the conformational transition. (a) The distribution of the autocorrelation function for the closed state in the Cartesian coordinates system. The single-exponential function is applied to fit the data for comparing. (b) The distribution of autocorrelation function for the closed state in the semilogarithmic coordinates system, in which the y value is logarithmic and the x value remains linear. A multiexponential function is applied to fit the data. (c) The distribution of the autocorrelation function for the open state in the Cartesian coordinates system with the single-exponential function fitting. (d) The distribution of the autocorrelation function for the open state in the semilogarithmic coordinates system, in which the y value is logarithmic and the x value remains linear. A multiexponential function is applied to fit the data.
maps. Each point in the contact maps represents one contact between two residues. A point with a red color means the that corresponding contact is near its closed native contact, and a point with a blue color implies the corresponding contact is close to its open native contact. The left column represent the typical corresponding structures. The distances of three pairs of residues (Asp10Lys115, Phe50Thr118, and Gly68Asp157) are shown in these figures to illustrate the meaning of points in the contact maps. Asp10, Lys115, Phe50, Thr118, Gly68, and Asp157 are labeled by D10, K115, F50, T118, G68, and D157, respectively. In the contact maps, each plot has three points enclosed by rectangles. They correspond to the three contacts (D10K115, F50T118, and G68D157). The arrows from the structures to the contact maps reflect the corresponding relationships between the distances and the points enclosed by rectangles. From the open state to the closed state, the distances of two contacts tend to decrease, and the color of the corresponding points in the contact maps change from blue to red. To emphasize the structural evolution of the open state, transition state, and closed state, we show those three contacts maps and their related structures in Figure 8. From Figure 811, we can see that the interaction between two groups of residues, which are located in region 1075 and region 115160, changes significantly when conformational transition occurs. Although the distances between these two groups of residues have some fluctuations during the conformational transition, group 1075 and group 115160 tend to get close to each other and form “forceps” of the protein to “clamp” the ligand when it comes in. Importantly, there are some
interesting residues in these two groups. Our simulations show that three residues, Asp157, Thr70, and Gly68, which are important for stabilizing the R-amino group of the ligand,9 move close when conformational transition occurs. The simulations also show that Asp157 moves faster to Thr70 than Gly68. Additionally, the residues Thr70, Arg75, and Gly119 interacting with the R-carboxyl group of the ligands also migrate together. Arg75 and Thr70 approach Gly119 rapidly. At the transition state, they have formed contacts which are very similar to those of the closed state (data not shown). An interesting characteristic of binding is the “doorkeeper”9 formed by residues Asp10 and Lys115, which locks the ligand within the binding pocket. The simulations also capture this feature. Additionally, Asp10 may forms a contact with Asn138 as well in the conformational transition. Another interesting feature of the contact formation is that from the open state to the closed state, most of the important contacts form slowly before the transition state, and after the transition state, the contacts form more quickly. We analyzed the φ values of each residues of the protein in the transition state, as shown in Figure 12. To be clear, we only exhibited the residues with high φ values. In the experiment, φ values provide an approach to quantify the strength of native interactions in the transition state.18,19 By calculating the φ value of a specific residue, one may know its importance in the protein dynamic processes, such as protein folding, conformational transition, and so on. However, it should be noted that φ values obtained from simulations of the simplified model may not correlate well with the experimental values. In our work, the 4153
dx.doi.org/10.1021/jp110845u |J. Phys. Chem. B 2011, 115, 4147–4159
The Journal of Physical Chemistry B
ARTICLE
Figure 8. Residues contact maps and the corresponding structures on the pathway from the open state to the closed state. In (a), (c), and (e), the distances of the residues (D10K115) are 23, 20, and 12.5 Å, respectively, the distances of the residues (F50T118) are 17.1, 16, and 7.8 Å, respectively, and for the residues (G68D157), the distances are 12.1, 11.6, and 10.1 Å, respectively. The colors of the corresponding points in the contact maps are from blue to light blue. For comparison, the points under the diagonal in each contact map represent the contacts of the structure of the closed state.
following formula20 was used to calculate the simulated φ values φi ¼
Æni ætran Æni æopen Æni æclosed Æni æopen
ð4Þ
where Æniæ is the average value of the number of contacts for residue i, and tran, open, and closed subscripts represent the
transition state, open state, and closed state, respectively. We can see that some important residues have a high φ value at the transition state. Thr70, Arg75, and Gly119 play an important role in stabilizing the ligand, and their φ values are 0.74, 0.77, and 0.71, respectively. Gln183 and Tyr185 from the second hinge may also participate in the conformational transition. Their φ values are 0.80 and 0.79, respectively. Other residues with 4154
dx.doi.org/10.1021/jp110845u |J. Phys. Chem. B 2011, 115, 4147–4159
The Journal of Physical Chemistry B
ARTICLE
Figure 9. Residues contact map and corresponding structures of structures 2, 3, and 4. In (a), (c), and (e), the distances of the residues (D10K115) are 20.8, 22.8 and 21.5 Å, respectively, the distances of the residues (F50T118) are 17.8, 15.3, and 15.3 Å, respectively, and for the residues (G68D157), the distances are 9.6, 10.0, and 12.4 Å, respectively. The colors of the corresponding points in the contact maps are from blue to light blue.
high φ values, such as Thr72 (φ value of 0.88), Tyr86 (φ value of 0.86), Asp122 (φ value of 0.87), Leu162 (φ value of 0.75), and Lys166 (φ value of 0.80), can be observed in the simulations, implying that these residues may be involved in the conformational transition as well and play important roles. We also can see that although some of residues, like Asp157, Gly68, Phe13, Phe50, Asp10, Lys115, and so on, are
important for stabilizing the ligand, they may not have high φ values.
’ DISCUSSION AND CONCLUSION Recently, some other types of two-well models have been developed to study the conformational transition.1113 Kei-ichi Okazaki et al. considered the topological characteristics of the 4155
dx.doi.org/10.1021/jp110845u |J. Phys. Chem. B 2011, 115, 4147–4159
The Journal of Physical Chemistry B
ARTICLE
Figure 10. Residues contact map and corresponding structures of structures 5, 7, and 8. In this set of contact maps, the colors of the corresponding points change from light blue to light green and light red. In (a), (c), and (e), the distances of the residues (D10K115) are 17.1, 20.7, and 13.1 Å, respectively, the distances between residues F50 and T118 are 20.0, 15.9, and 11.1 Å, respectively; and for residues G68 and D157, the distances change from 12 to 10.7 Å.
energy landscape of proteins and created two independent structure-based potentials and connected them smoothly to make a double-well model. The same spirit can be found in the work of Best et al. Comparing these models, one feature of our model is to consider the microscopic situation of each contact pair. According to the known structures, we specified the form of interaction for each pair. Therefore, we call other types of models macroscopic two-well models. We call our models microscopic
two-well models. While both models give stable open and closed states reflected in free-energy profiles, the macroscopic model leads to the macroscopic two wells, and but it is hard to get prominent intermediate states. The microscopic models can lead to intermediate states between open and closed states. Both models represent some extreme situations. Macroscopic models more likely represent the cases with smooth underlying energy landscapes, while microscopic models more likely represent the 4156
dx.doi.org/10.1021/jp110845u |J. Phys. Chem. B 2011, 115, 4147–4159
The Journal of Physical Chemistry B
ARTICLE
Figure 11. Residues contact map and corresponding structures of structures 9 and 10. In this set of contact maps, the color of the corresponding points becomes deep red. In (a) and (c), the distances of the residues (D10K115) are 14.8 and 13.3 Å, respectively, the distances between residues F50 and T118 are 9.4 and 8.9 Å, respectively, and for the residues G68 and D157, the distances change from 10.8 to 9.1 Å.
cases with more rugged underlying landscapes. It is worth mentioning that other kinds of two-well model were also developed to study the problem of conformational switch. In the work of Jhih-Wei Chu et al.,21 the authors connected two elastic network models to generate a double-well network model in which the authors considered the detailed interactions between residues with connecting two-harmonic potentials, and they applied this model to study the two conformational systems, that is, the coil-to-helix transition of the DB loop in G-action and the open-to-closed transition of adenylate kinase. The freeenergy function generated by this model can identify the distinct minimum free-energy paths between two conformational states. The model in the work of Kei-ichi Okazaki et al. has been further developed to deal with the effect of the ligand and the binding mechanism. They considered four states, that is, the unbound-open state, bound-open state, unboundclosed state, and bound-closed state, of the protein. In our model, the interaction between the ligand and the protein is included implicitly in the potential energy function because the native contacts of the closed state of the protein include the ligand, which may induce the formation of the native contacts of the closed state. In the future, we may further consider more details of the binding mechanism. In this paper, we developed a structure-based two-well model17 to study the properties of the kinetics and statistical distributions for the conformational transition of GlnBP, which is
one of the small periplasmic binding proteins that carry small ligands form the periplasmic space into the cytoplasmic space. In the process of conformational transition, GlnBP exhibits two stable states, that is, the ligand-free open state and ligand-bound closed state, implying that the protein may go through two distinct local minima on the potential surface. We constructed the free-energy landscapes of the conformational transition with different temperatures and analyzed their topological characteristics. Two basins were observed on the free-energy landscapes with the temperatures under the melting point. One corresponds to the open state, and the other represents the closed state. In our simulation, there was no detectable intermediate state in the process of conformational transition. The topological properties of the free-energy landscape show that the protein prefers to stay in the closed state at low temperature and tends to dwell in the open state at high temperatures. We also studied the mean and distribution of the first passage time. Both the closed- and opendwelling times exhibit the Γ distribution. With different temperatures, the scale parameters λo and λc of the Γ distributions decrease or increase monotonously. The analysis of the autocorrelation coefficient shows that the conformational transition may be a multiexponent process, which indicates the complexity and hierarchical structure of the underlying energy landscape. Finally, the contact maps and φ values of each residue were determined to illustrate the structural evolution and the important residues in the conformational transition. Some residues that 4157
dx.doi.org/10.1021/jp110845u |J. Phys. Chem. B 2011, 115, 4147–4159
The Journal of Physical Chemistry B
ARTICLE
Figure 12. (a) φ values of residues with a high value. (b) A typical structure in the transition state. The residues with high φ values in the large domain (Thr70, Thr72, Arg75), the small domain (Gly119, Asp122, Lys166), and the hinge region (Gln183, Try185, Tyr86) are marked in red.
are critical for binding and stabilizing the ligand show high φ values and significant transient between the open state and the closed state.
’ AUTHOR INFORMATION Corresponding Author
*E-mail:
[email protected].
’ ACKNOWLEDGMENT The authors thank the National Science Foundation Career Award and National Institute of Health for support. ’ REFERENCES (1) Royer, C. A. Chem. Rev. 2006, 106, 1769–1784. (2) Lillo, M. P.; Beechem, J. M.; Szpikowska, B. K.; Sherman, M. A.; Mas, M. T. Biochemistry 1997, 36, 11261.
(3) Levy, Y.; Wolynes, P. G.; Onuchic, J. N. Proc. Natl. Acad. Sci. U.S.A. 2004, 101, 511–516. (4) Shoemaker, B. A.; Portman, J. J.; Wolynes, P. G. Proc. Natl. Acad. Sci. U.S.A. 2000, 97, 8868–8873. (5) Papoian, G. A.; Wolynes, P. G. Biopolymers 2003, 68, 333–349. (6) Wang, J.; Verkhivker, G. M. Phys. Rev. Lett. 2003, 90, 188101. (7) Go, N. Annu. Rev. Biophys. Bioeng. 1983, 12, 183–210. (8) Hsiao, C. D.; Sun, Y. J.; Wang, B. C. J. Mol. Biol. 1996, 262, 225–242. (9) Sun, Y. J.; Rose, J.; Wang, B. C.; Hsiao, C. D. J. Mol. Biol. 1998, 278, 219–229. (10) Cheung, M. S.; Garcia, A. E.; Onuchic, J. N. Proc. Natl. Acad. Sci. U.S.A. 2002, 99, 685–690. (11) Best, R. B.; Chen, Y. G.; Hummer, G. Structure 2005, 13, 1755–173. (12) Okazaki, K. I.; Koga, N.; Takada, S.; Onuchic, J. N.; Wolynes, P. G. Proc. Natl. Acad. Sci. U.S.A. 2006, 103, 11844–11849. (13) Okazaki, K. I.; Takada, S. Proc. Natl. Acad. Sci. U.S.A. 2008, 105, 11182–11187. 4158
dx.doi.org/10.1021/jp110845u |J. Phys. Chem. B 2011, 115, 4147–4159
The Journal of Physical Chemistry B
ARTICLE
(14) Sobolev, V.; Wade, R. C.; Vnend, G.; Edelman, M. Proteins: Struct., Funct., Genet. 1996, 25, 120–129. (15) Edman, L.; Mets, U.; Rigler, R. Proc. Natl. Acad. Sci. U.S.A. 1996, 93, 6710–6715. (16) Palmer, R. G.; Stein, D. L.; Abrahams, E.; Anderson, P. W. Phys. Rev. Lett. 1984, 54, 958–961. (17) Lu, Q.; Wang, J. J. Am. Chem. Soc. 2008, 130, 4772–4783. (18) Fersht, A. R.; Leatherbarrow, R. J.; Wells, T. N. C. Biochemistry 1987, 26, 6030–6038. (19) Matouschek, A.; Fersht, A. R. Methods Enzymol. 1991, 202, 81–112. (20) Ejtehada, M. R.; Avall, S. P.; Plotkin, S. S. Proc. Natl. Acad. Sci. U.S.A. 2004, 101, 15088–15093. (21) Chu, J. W.; Voth, G. A Biophys. J. 2007, 93, 3860–3871.
4159
dx.doi.org/10.1021/jp110845u |J. Phys. Chem. B 2011, 115, 4147–4159