3060
J. Phys. Chem. B 2010, 114, 3060–3069
Crucial Roles of the Subnanosecond Local Dynamics of the Flap Tips in the Global Conformational Changes of HIV-1 Protease Dechang Li,† Baohua Ji,*,‡ Kehchih Hwang,† and Yonggang Huang*,§ Department of Engineering Mechanics, School of Aerospace, Tsinghua UniVersity, Beijing 100084, China, Biomechanics and Biomaterials Laboratory, Department of Applied Mechanics, Beijing Institute of Technology, Beijing 100081, China, and Department of CiVil and EnVironmental Engineering, Northwestern UniVersity, EVanston, Illinois 60208 ReceiVed: January 20, 2010
To understand the underlying mechanisms of the open and closed conformational change of HIV-1 protease (HIV-1 PR) at multiple time scales, we performed serial fully unrestrained, extremely long time molecular dynamics simulations with an explicit solvent model. Spontaneous semiopen to closed conformational transition and inhibitor-collision-induced opening of the flaps were simulated in a real time scale. We found that the rapid, local subnanosecond fluctuations of the flap tips might be the mechanisms triggering the global open and close conformational transitions at the 100-ns time scale. The subnanosecond fluctuation is induced by the Φ-Ψ rotations of the residues at the flap tips, mainly Ψ of Gly49 and Φ of Ile50, disturbing the interactions between the two tips and then their stability. We further showed that the water molecule W301 is helpful for the stability of the PR-inhibitor complex by acting as a collision buffer for the dynamic interaction between flap tips and the inhibitor. These results might help gain a better insight into the dynamics of HIV-1 PR, especially the local dynamics of the flap tips, which may provide important guidelines for design of novel potent inhibitors. 1. Introduction Due to its central role in processing viral maturation, HIV-1 protease (HIV-1 PR) continues to be one of the primary targets of anti-AIDS therapy.1-4 Current inhibitors have limited longterm use because of the selection of inhibitor-resistant variants of the protease.5-7 HIV-1 PR is a dimeric aspartic protease in which each monomer contains 99 residues and the active site is caved by two flexible β-hairpin flaps (i.e., residues 43-58). Crystallographic and NMR studies have showed that the flaps exist in a large ensemble of conformations, mainly categorized by semiopen, opened, and closed states, which are closely related to the gated ligand binding processes, as shown by Figure 1. Furthermore, upon binding with the protease, the inhibitors formed hydrogen bonds either directly or mediated by a water molecule, W301, with the flap tips,8-11 stabilizing the conformation of the flaps and then the bound protease complex. Therefore, the understanding of the dynamics of the flaps is critically important for inhibitor designs. NMR studies have indicated that the entire flaps move at a 100 µs time scale, but the flap tips undergo fast fluctuation on a subnanosecond time scale for free HIV-1 PR.12-14 However, the dynamics of the subnanosecond fluctuation (local) of the flap tips and its relation with the open and close conformational change (global) of flaps at large time scale has not been understood. In addition to the experimental studies, lots of efforts in molecular dynamics (MD) simulations have been paid to HIV-1 PR to reveal the underlying mechanisms of its dynamics motions. Scott and Schifer15 simulated the spontaneous opening * Corresponding authors. (B.J.) Phone: +86-10-68918309. E-mail:
[email protected]. (Y.H.) Phone: (847) 467-3165. E-mail: y-huang@ northwestern.edu. † Tsinghua University. ‡ Beijing Institute of Technology. § Northwestern University.
of the flaps of free protease with a 10 ns MD simulation. Their results suggested that the curling of flap tips acts as a trigger for the wide opening of flaps. However, Meagher and Carlson16 doubted their simulations by claiming that the fast opening events may be caused by insufficient pre-equilibrium. Long time scale simulations with sufficient pre-equilibrium steps are necessary for reliable simulation results. Later, To´th and Borics17 reported a single opening event based on their restrained MD simulations with implicit aqueous solvent in which the bond lengths and bond angles were set rigid. They showed that there was a network of weak polar interactions between the flap tips, which was proposed to be highly sensitive to the flap fluctuations.17,18 It was shown that the flaps can be opened at about 100 ps time scale by a small external impulse near the flap region.19 The first report about both open and reclosing events of the flaps in one simulation were given by Hornak et al.20 using implicit aqueous solvent without any restraints of the system. However, the use of artificially low solvent viscosity precludes the direct comparison of their simulation with the NMR experiments because of a large difference in the time scale between their simulations and the NMR experiments. In addition, the implicit solvent model lacks details of the water molecules; therefore, it cannot consider the effect of a special single water molecule, for example, W301, which was found in most of the inhibitor bound complexes. Recently, Ding et al.21 performed a MD simulation with an explicit solvent model to obtain the distance distribution of the labels linked at residues 55/55′. They showed that the simulations with the explicit solvent model can be directly compared with the experiments. At the same time, acceleration methods were proposed to simulate the large conformational change of protease. Hamilberg et al.22 and Wiley et al.23 developed accelerated MD simulation methods with an explicit solvent model to investigate the dynamic motions of the flaps. These accelerated simulations
10.1021/jp1005549 2010 American Chemical Society Published on Web 02/09/2010
The Roles of Local Dynamics at the Flap Tips
J. Phys. Chem. B, Vol. 114, No. 8, 2010 3061
Figure 1. Cartoon draws of three distinct conformations of HIV-1 protease: (a) semiopen, (b) closed, and (c) opened. The flap tips are highlighted in red; the rest of the flaps are in yellow. The residues Asp located at the active site are displayed in VDW spheres. There is an inhibitor XK263 filled in the active cavity of the closed conformation, represented by green rods. Top views illustrate the distance between and handedness of the flaps in the three conformations.
were achieved by modifying either the force field or the algorithms to enhance the fluctuations of the flap tips to induce the opening events. These methods were supposed to bring artificial motions to the system. On the other hand, McCammon and co-workers24,25 developed a coarse-grained model of HIV-1 PR, by which multiple opening and closing events were simulated in one simulation. Their model was successfully used for simulating the ligand binding and product release processes.26-29 However, the current coarse-grained models lacked atomic details about the crucial local structure and interaction, such as the special hydrogen bond networks between protease and ligands, which should be important for inhibitor designs. Although previous simulations have provided a variety of thought-provoking views of the dynamics of the protease, their results cannot straightforwardly be compared with the experimental measurements due to the limitations of their approaches. In addition, these simulations were mostly focused on the entire motions of the flaps, but the subnanosecond time scale motions of the tips were rarely reported. Therefore, long time scale, allatom MD simulations with an explicit solvent model are highly desirable to investigate the dynamics of HIV-1 PR in multiple time scales, that is, from subnanosecond to 100 ns. In this work, we applied serial fully unrestrained, extremely long MD simulations (total ∼1.62 µs, see Table S1 of Supporting Information) with an explicit solvent model to investigate the underlying mechanisms of the dynamics of HIV-1 PR. We simulated the conformational transitions among the semiopen, opened, and closed states. We found that the subnanosecond local motion of the flap tips is crucial to the global conformational change of the protease. This subnanosecond fluctuation in our simulations was suggested to correspond to the subnanosecond motion observed in the NMR experiments. 2. Methods 2.1. MD Simulations with Explicit Solvent. The MD simulations were performed using the Gromacs package30,31 with the AMBER force field of ffamber99,32 in which the all-atom force field parameters of inhibitors were obtained by the ANTECHAMBER module and GAFF33,34 with AM1-BCC35 charges in the AMBER package.36 The system was solvated in a 90 × 80 × 80 Å3 TIP3P37 water box, with about 15 000 water molecules. Appropriate chlorine ions were added to neutralize the system. The particle mesh Ewald (PME)38 method was used to calculate the long-range electrostatic interactions. To avoid insufficient pre-equilibration,16 the systems were minimized first using a steepest descent
algorithm of 10 000 steps. Then, the system was gradually heated to 300 K in 200 ps while positional restraints were used and the restraint force constants were gradually decreased from 2.39 to 0 kcal/(mol · Å2) in a few stages. All production simulations were fully unrestrained at 300 K with a pressure of 1 bar with the Berendsen algorithm.39 The SHAKE algorithm40 was applied to constrain the bonds with H-atoms. The time step of the simulations was 2.0 fs. The cutoff of the nonbonded interactions was set to 10 Å. The nonbonded pairs were updated in every 10 steps. All graphics and visualization analysis were processed using the VMD program.41 2.2. Preparation of Models. The structures for protease were retrieved from the Protein Data Bank with PDB codes: 1HHP42 for free protease of semiopen conformation (Figure 1a), 1HVR8 for the free protease of closed conformation by deleting the ligand XK263 (Figure 1b), 1HVR8 for the XK263-bound complex, 1MES9 for the DMP323-bound complex, 1HXB10 for the Saquinavir (SQV)-bound complex, and 1HXW11 for the Ritonavir (RTV)-bound complex. The chemical structures of the inhibitors can be found in Figure S1 of the Supporting Information. A widely opened structure (Figure 1c) of the protease was retrieved from our full-atom MD simulation induced by collision with inhibitor XK263 from outside. The side chains of the atomistic structure of proteases except the one from IHHP were modified using the WHAT IF web interface43 to match the sequence of the PDB of 1HHP. The catalytic Asp side chains of the free protease without bound ligands were modeled in the nonprotonated state according to the protonation in solution at pH 7, whereas for the bound complex, the catalytic Asp side chains were protonated according to the experiments and calculations; that is, both side chains of Asp25/Asp25′ were protonated for DMP323 and XK263 bound complex,44 and one Asp25/Asp25′ was protonated for RTV and SQV bound complex.45-47 3. Results and Discussion 3.1. Global Conformational Change of Flaps. To study the dynamics of the conformational change of HIV-1 PR, we used different initial structures for the simulations; that is, semiopen, closed, and opened conformations, with or without interaction with an inhibitor. The preparation of models can be found in the Methods section. The conformational changes of the protease were monitored by the root-mean-square deviation (rmsd) of CR atoms of the flaps (residues 43-58) with respect to two reference structures: the semiopen and closed structures. Figure 2a shows the simulation starting from the semiopen structure,
3062
J. Phys. Chem. B, Vol. 114, No. 8, 2010
Li et al.
Figure 2. Root mean square deviation of flap CR atoms in simulations starting from (a) semiopen and (b) closed conformations, as well as (c) the simulation from the semiopen conformation, in which the protease interacts with an inhibitor, XK263, outside the active site. The decrease in the rmsd of the red curve indicates the structure is transforming to the closed structure, and the decrease in the rmsd of the black curve indicates the structure is transforming to the semiopen structure. The arrows highlight the events of conformational transitions.
Figure 2b shows the simulation starting from the closed structure, and Figure 2c shows the simulation from the semiopen structure interacting with inhibitor XK263. The decrease in the rmsd of the red curve indicates the structure is transforming to the closed structure, and the decrease in rmsd of the black curve indicates the structure is transforming to the semiopen structure. We can see that the two rmsd curves in Figure 2a meet at t ∼ 125 ns and cross over each other after ∼10 ns. These results indicate that the structure is transforming from the semiopen to the closed conformation. At the beginning of the transition, one flap may go up above the other, as shown by Figure 3a at t ) 101.88 ns, which subsequently causes the handedness switch of the flaps, as shown by Figure 3a at t ) 102.50 and 128.33 ns. The simulations also show that the transition from closed to semiopen structures is happening at t ∼ 200 ns, indicated by the second meeting of the two rmsd curves of Figure 2a as well as t ∼ 80 ns of Figure 2b. However, the flaps then returned quickly and resumed the handedness of the closed conformation
in a short time. More details of the conformational snapshots during the simulation can be found in Figure S2 of the Supporting Information. Figure 2c shows a semiopen structure can transform to an open conformation by the collision of an inhibitor. The value of the two rmsd curves increases to as large as 10 Å, indicating that the structure deviates from both semiopen and closed conformations. Figure 3b shows that the collision of the inhibitor can disturb the interaction between the two flap tips. Once the interaction between the flap tips is weakened and the distance between the tips is enlarged beyond a threshold value, the semiopen-to-open conformational transition occurred, as shown by Figure 3b at t ) 50 ns. The above results show that the protease can perform conformational changes that are closely related to the dynamics of the flaps. The intrinsic fluctuation of the flaps and the collision of the inhibitor would influence the interaction between the flap tips (Figure 4). NMR studies have showed that the flap tips undergo subnanosecond time scale dynamics.13,14 The biological
The Roles of Local Dynamics at the Flap Tips
J. Phys. Chem. B, Vol. 114, No. 8, 2010 3063
Figure 3. (a) Snapshots of the conformational changes of HIV-1 PR during the simulation starting from the semiopen conformation. (b) Snapshots of the conformational change of HIV-1 PR interacting with inhibitor XK263 outside the active site. The flaps are blue, the tips are red, and the flap elbow is pink. The residues Asp located at the active site are displayed in VDW spheres. The inhibitor XK263 yellow in the VDW spheres.
Figure 4. Illustration of interacting dynamics between the flap tips of free protease (a) without interaction of inhibitor and (b) with interaction of the inhibitor outside the active cavity. The inhibitor XK263 is yellow in the VDW spheres. The numbers in the graph indicate the distances between the two tips.
implications of the local conformation change at the flap tips at the nanosecond scale will be studied in the following. 3.2. Local Dynamics of the Flap Tips. The experiments suggested that the flap tips undergo subnanosecond time scale fluctuation in the free protease.12-14 To determine the modes of the subnanosecond motions, the hydrogen bond length at the flaps was calculated for analyzing the dynamics of the flap tips. Figure 5a and b shows the evolution of the intraflap hydrogen bond length of the two flap tips during the simulation. As we can see, the intraflap hydrogen bonds at the flap tips between Gly49 and Gly52 are less stable as compared with the other hydrogen bonds (Figure 5a and b, and Figure S3 in the Supporting Information), with the bond length frequently increased to as large as 4-6 Å, indicating the rupture of the bonds. This result is consistent with the cross-strand amide nuclear Overhauser effect data.13 Figure 5c and d shows that the change of the bond length between Gly52 and Gly49 is caused by rotating dynamics of the dihedral angle Ψ of Gly49 at the tip. Along with the Ψ rotation of Gly49, the dihedral angle Φ of Ile50 will change simultaneously. These rotations will also change the interaction between the two tips, which may further influence the stability of intraflap hydrogen bonds of the tips. For example, the rotations cause the rupture of intraflap hydrogen bonds as well as the formation of one
interflap hydrogen bond between the two tips (see Figure 5e and Figure S4 of the Supporting Information) that restrained the Φ-Ψ rotations so that the intraflap hydrogen bonds could not reform during 60-100 ns (see Figure 5a and b). Once the interflap hydrogen bond between the tips is ruptured, the intraflap hydrogen bonds are reformed by the Φ-Ψ rotations, triggering the conformational transitions of the flaps, that is, the switch of the handedness (see Figure 2a). We then calculated the persistent and broken times of the intraflap hydrogen bonds at the flap tips, as shown in Table 1. The persistent times of the intraflap hydrogen bonds at the tips were at nanosecond time scale, and the broken times were at subnanosecond time scale. These results indicate that the subnanosecond time scale motions of flap tips observed in the experiments correspond to the breaking of the intraflap hydrogen bonds, which is induced by the Φ-Ψ rotations of the residues. In comparison, the persistent times of the intraflap hydrogen bonds of opened conformation are much shorter than those of the semiopen and closed ones (i.e. ∼1 ns), which might be caused by the absence of the interaction between the two flap tips due to their large separation. We further studied the change of the angles formed by the triplet CR atoms of the tips, that is, residues Gly48-Gly49-Ile50 and Gly49-Ile50-Gly51. The distributions (probabilities) of the angles in the simulations started from semiopen and opened conformations were calculated. Figure 6 shows that the two angles had bimodal distributions in the semiopen and opened conformations, which are consistent with previous simulations.48 The insets in Figure 6 show that the bimodal distributions are closely related to the Φ-Ψ rotation dynamics of the tips. It is noted that one of the peaks is much lower than the other in the semiopen conformation, but the two peaks are similar in height in the opened one. The lower peak corresponds to the state with a broken intraflap hydrogen bond within flap tip; the higher peak, to the state with an intact intraflap hydrogen bond. These observations are consistent with our simulations that the persistent times of the hydrogen bonds at the tips are much larger than the broken times in the semiopen conformation, and the two times are roughly the same in the opened conformation. This might be the physical basis for the quadratic double-welltype angle potential in the coarse-grained model developed by McCammon and co-workers.24,25 Our simulations also show that the change of the angles formed by the triplet CR atoms can
3064
J. Phys. Chem. B, Vol. 114, No. 8, 2010
Li et al.
Figure 5. The evolution of the intraflap hydrogen bond length of the tips of (a) flap A and (b) flap B during the simulation starting from the semiopen conformation. The red curve is the average values with a running average time window of 100 ps. (c) The hydrogen bonds are presented between Gly49 and Gly52. Dihedral angles Φ, Ψ, and ω are illustrated by arrows. (d) The hydrogen bonds of the tips disrupted by the Φ-Ψ rotations. (e) Illustration of hydrogen bonds formed between the two tips because of the Φ-Ψ rotations.
TABLE 1: The Persistent and Broken Times of the Intraflap Hydrogen Bonds at Flap Tipsa semiopen conformation 52O-49N 52N-49O 52′O-49′N 52′N-49′O
closed conformation
opened conformation
persistent time (ps)
broken time (ps)
persistent time (ps)
broken time (ps)
persistent time (ps)
broken time (ps)
7233 ( 9073 2680 ( 4077 6768 ( 10182 6318 ( 7444
734 ( 1597 240 ( 492 1373 ( 2648 653 ( 846
10018 ( 15586 4144 ( 7303 10363 ( 13499 10252 ( 12024
536 ( 660 344 ( 540 791 ( 998 859 ( 1107
591 ( 1273 1739 ( 1944 500 ( 828 948 ( 956
520 ( 866 985 ( 789 691 ( 850 1615 ( 1748
a The persistent time is defined as the average time that the hydrogen bond is intact; the broken time is defined as the average time that the hydrogen bond is broken.
cause the curling of the flap tips. Scott and Schiffer15 suggested that the curling was related to the flaps opening, according to their study. Wiley et al.23 reported an accelerated MD simulation of the HIV-1 PR system in which the residues located at the flap tips (i.e., 49-51) were “heated” to temperatures as high as ∼1000 K. The high temperature enhanced the fluctuation motions at the flap tips, which caused the opening of the flaps and transitions between the semiopen and closed conformations. These studies suggested that the fluctuation
motions of the flap tips are closely related to the conformational transitions. Hamelberg and McCammon22 reported another accelerated MD simulation of the protease in which the potential surface of system was modified to enhance the escape rate from one potential well to another. They suggested that the local dynamics of the flap tips is the underlying mechanism that causes the opening of the flaps. In their simulations, they observed frequent trans-cis isomerization of the Gly-Gly ω-bonds and proposed that the isomerization might be related to the opening events.
The Roles of Local Dynamics at the Flap Tips
Figure 6. The distributions of the angles formed by the continuous triplet CR atoms of the flap tips (a) residue Gly48-Gly49-Ile50 and (b) residue Gly49-Ile50-Gly51, from the simulations starting from both semiopen and opened conformations. The insets illustrate how the Φ-Ψ rotations change the angle of the triplet CR atoms.
However, we did not observe any such isomerization in our simulations. Instead, we saw frequent dihedral angle (Φ-Ψ) rotations. Figure 7 shows Φ-Ψ free-energy profiles of the residues on the tips from the simulations of both semiopen and opened conformations. There are multiminima of the free energy profiles for most of the residues, especially for Gly49 and Ile50. The energy barrier of Φ-Ψ rotations is about 6kBT, according to our calculations, whereas that of the trans-cis isomerization of the Gly-Gly ω-bonds is about 12kBT.22 Therefore, the Φ-Ψ rotations should happen more easily and have a higher probability than the transition of the ω-bonds. The above analyses suggest that the subnanosecond dynamics of the flaps are induced by the Φ-Ψ rotations of the residues at the flap tips. The Φ-Ψ rotations disturb the interaction between flap tips and induce curling of the tips, which consequently causes the conformational changes of the protease. 3.3. Comparison of Simulations with NMR Experiments on the Local Dynamics of Flap Tips. To compare the characteristics of the local dynamics of HIV-1 protease inferred from MD simulation with the experimental data,13,49,50 we calculated the generalized order parameters S2 of the backbone amide, a measure of the spatial restriction of the internal motions, for
J. Phys. Chem. B, Vol. 114, No. 8, 2010 3065 free and ligand-bound protease and the effective internal correlation times τe, a measure of the rate of the motions, from the so-called model-free approach.50 τe can be the measure of the time scale of the local dynamics at the flap tips. The calculated results of S2 agree with the experiment data,13,49 as shown in Figure 8a and b. Figure 8a shows that residues 15-17, 37-41 (flap elbow), 49-53 (flap tips), 68-71, and Thr80 in the free protease are more flexible. In the flap tips region, S2 ∼ 0.5 at residue Ile50, indicating that its NH bond reorients in a range of as much as 38° (using the cone model to interpret the physics of S2). In contrast, the flexibility of the flap tip region is suppressed by the inhibitor in the bound protease, as shown in Figure 8b. It also shows that the flap tips in the bound complex with inhibitor RTV are more stable than those with inhibitor DMP323. The effective internal correlation times, τe, calculated from oursimulationswerealsoconsistentwiththeNMRexperiments.13,49,51 Figure 8c shows the results for the free protease, indicating that the effective correlation times, τe, of the motions at the region of flap tips are at the subnanosecond time scale. In contrast, the value of τe at the flap tips of the inhibitor-bound complexes are much shorter, being a few picoseconds, as shown in Figure 8d. These results suggest that our simulations have captured the main dynamics of local motion of the flaps found in experiments in a real time scale. The large differences in the dynamics of the flap tips between the free protease and the bound complex suggest the important biological function of the flap tips. 3.4. Effect of Water Molecule W301 on the Dynamics of Flap Tips. In the bound complex, the hydrogen bonds either form directly (e.g., DMP323, XK263) or are mediated by a water molecule (e.g., RTV, SQV) between the flap tips and the inhibitor. In the DMP323/XK263 bound complex, residue Ile50 at the flap tips forms hydrogen bonds with the urea oxygen of the inhibitor. When the Φ-Ψ rotations of the flap tips occur, they can disrupt the hydrogen bonds between the tips and urea oxygen and then cause instability of the flap tips, as shown in Figure 9a. In contrast, in the RTV/SQV bound complex, the hydrogen bond between flap tips and the inhibitor is mediated by a water molecule (W301). The Ile50-W301-inhibitor serial hydrogen bonds are much more flexible. When the Φ-Ψ rotations happen, the water molecule may arrange its position and pose to fit the subtle conformational changes of the tips (see Figure 9b). Therefore, the water molecule can be helpful for the stability of the hydrogen bond networks between flap tips and inhibitor. To confirm this hypothesis, we calculated the average hydrogen bond lengths and angles between the flap tips and the inhibitors for the XK263/DMP323 bound complex as well as those between the flap tips and W301 and those between W301 and inhibitor for the RTV/SQV bound complex, as shown in Table 2. The results indicate that the hydrogen bond lifetime between the flap tips and W301 are longer than that between the flap tips and DMP323/XK263 (see Table 2), which is consistent with the quantum mechanical calculation results.52 Furthermore, Figure 10a shows that the distributions of the hydrogen bond length between the flap tips and DMP323/ XK263 display a wider range than those between the tips and W301, which indicates that the hydrogen bond between the flap tips and W301 are more stable. It is noteworthy that the distribution of the hydrogendonor-acceptor angle between W301 and SQV exhibits a bimodal shape without hydrogen bonds’ rupture, as shown in Figure 10b. To understand this, the position and post of the
3066
J. Phys. Chem. B, Vol. 114, No. 8, 2010
Li et al.
Figure 7. Φ-Ψ free-energy profiles of the residues at the tips from the simulations starting from (a) semiopen and (b) opened conformations. The free energy profiles were obtained by the Boltzmann inversion method. The unit of the color bar is kBT.
W301 are illustrated as the inset in Figure 10b. As we can see, W301 forms two hydrogen bonds with the CO groups at the
P1′ and P2 sites of SQV, respectively. W301 can rotate around the C2 symmetry axis of the protease to switch the position of
The Roles of Local Dynamics at the Flap Tips
J. Phys. Chem. B, Vol. 114, No. 8, 2010 3067 also supported by the experimental observations that the RTV and SQV inhibitors have a much lower dissociation rate (∼105 fold) than that of the cyclic urea analog inhibitors.53 In fact, most of the inhibitors for HIV-1 protease have a water molecule with the flap tips in the crystal structures. 4. Conclusions
Figure 8. The generalized order parameters S2 for (a) free and (b) inhibitor-bound protease and the effective internal correlation times, τe, for (c) free and (d) inhibitor-bound protease. The experimental data were obtained from previous published works.13,49,51 S2 and τe were calculated by averaging the values of the two monomers, with the error bars reflecting the differences between the monomers.
its two hydrogen atoms between the hydrogen bonds at the P1′ and P2 sites, while the integrity of the hydrogen bond networks are still maintained between the flap tips and W301. These observations support our hypothesis that W301 serves as a collision buffer in the hydrogen bond networks. Our results are
In this study, we performed serial fully unrestrained, longtime, all-atom MD simulations with an explicit solvent model on HIV-1 PR systems starting from different conformations; that is, semiopen, closed, and opened states of the protease. The local dynamics of the flap tips were analyzed during the global conformational change of the protease from semiopen to closed and semiopen to open states. Our simulations showed that the subnanosecond fluctuation induced by the Φ-Ψ rotations of the residues at the tips can trigger the open and closed states of the flaps. The model-free analysis demonstrated that the dynamic characteristics of the subnanosecond motions at the flap tips obtained from our simulations were in good agreement with the results from the NMR experiments.12-14 The underlying physics in the relationship between the subnanosecond local dynamics at the flap tips and the global dynamics of the flaps can be understood as follows. The Φ-Ψ rotations of the residues at the tips can disturb the weak interactions between the two tips and then their stability, which consequently causes the global conformational transition of the flaps. We calculated the free energy profiles of the Φ-Ψ rotations of the flap tip residues. The results showed that the energy barriers of the Φ-Ψ rotations are not higher than 6kBT, which is much lower than that of the Gly-Gly ω-bond transitions.22 Therefore, the Φ-Ψ rotations of flap tips are more likely to be the trigger of the global conformational transitions of flaps instead of the Gly-Gly ω-bonds transitions. There may be other factors that can influence the conformational transition of the flaps, such as the interaction of the flaps with the rest of the protease. For example, we found that there was one heteromonomer hydrophobic cluster formed by the tip (Ile50) of one monomer with the P1-loop (residues 79-81) of the other one, which made the two tips slightly separate (see Figure S5 of Supporting Information), indicated by a short-lived jumping at about 247 ns in Figure 2a, which is consistent with previous studies.54,55 We further found that the hydrogen bonds formed between flap tips and W301 are more stable than those formed directly between flap tips and the inhibitors. This is because water
Figure 9. The hydrogen bond networks between the inhibitor and flap tips of (a) XK263 and (b) SQV bound complex. The numbers illustrate the length of the hydrogen bonds.
3068
J. Phys. Chem. B, Vol. 114, No. 8, 2010
Li et al.
TABLE 2: The Properties of Hydrogen Bonds between Flap Tips and Inhibitors, Those between the Tips and W301, and Those between W301 and the Inhibitor crystal structure bound complex XK263 DMP323 SQV RTV a
flap tips-XK263 flap tips-DMP323 flap tips-W301 W301-SQV flap tips-W301 W301-RTV
dimulations
length (Å)
angle (°)
length (Å)
angle (°)
H-bond energy (kcal/mol)
H-bond lifetimea (ps)
3.21 ( 0.05 3.27 ( 0.05 3.07 ( 0.23 2.87 ( 0.04 2.97 ( 0.13 2.68 ( 0.04
17.26 ( 5.93 7.07 ( 5.19 n.a. n.a. 16.46 ( 13.65 10.09 ( 0.51
3.26 ( 0.40 3.30 ( 0.43 3.03 ( 0.20 2.90 ( 0.19 3.06 ( 0.18 2.74 ( 0.13
23.77 ( 23.25 25.59 ( 19.92 12.78 ( 6.96 12.76 ( 7.36 15.76 ( 9.36 10.26 ( 5.35
-1.65 ( 1.14 -1.51 ( 1.01 -2.34 ( 0.40 -8.20 ( 0.90 -2.18 ( 0.41 -9.32 ( 0.94
749.71 546.22 1684.55 523.17 2400.80 514.18
The H-bond lifetime was estimated according to Luzar and Chandler.56
structures of inhibitors used in the simulations; chronological order views of the procedure of flap handedness reversal in the simulation starting from semiopen conformation; plots of lengths of intraflap hydrogen bonds of the β-hairpin flaps, except the tips during the simulation from semiopen conformation; a figure of average numbers of hydrogen bonds formed between the two flap tips during the simulation starting from semiopen and closed conformations; a figure of the top view of the hydrophobic cluster formed by the tip (Ile50) with the P1 loop (residues 79-81) of the other monomer at 247.40 ns of the simulation starting from semiopen conformation.This material is available free of charge via the Internet at http://pubs.acs.org. References and Notes
Figure 10. Distributions of (a) hydrogen bond lengths and (b) hydrogen-donor-acceptor angles between the flap tips and inhibitor, and those between the flap tips and W301, as well as those between W301 and the inhibitor. The insets in part (b) illustrate the switch of the H atoms of W301 without breaking the hydrogen bonds between the W301 and SQV.
molecule W301 can maintain the integrity of the hydrogen bond networks between the tips and inhibitor through its flexible movement and rotations. We therefore suggest that the hydrogen bonds mediated by W301 between flap tips and inhibitors are in favor of the stability of the bound complex. Acknowledgment. This research was supported by the National Natural Science Foundation of China through Grants nos. 10628205, 10732050, 10872115, and the National Basic Research Program of China through Grants nos. 2004CB619304 and 2007CB936803. Supporting Information Available: A table for the length of simulation time of various systems; a figure of the chemical
(1) Kohl, N. E.; Emini, E. A.; Schleif, W. A.; Davis, L. J.; Heimbach, J. C.; Dixon, R. A. F.; Scolnick, E. M.; Sigal, I. S. Proc. Natl. Acad. Sci. U.S.A. 1988, 85, 4686. (2) Katz, R. A.; Skalka, A. M. Annu. ReV. Biochem. 1994, 63, 133. (3) Wlodawer, A.; Vondrasek, J. Annu. ReV. Biophys. Biomol. Struct. 1998, 27, 249. (4) Hornak, V.; Simmerling, C. Drug DiscoVery Today 2007, 12, 132. (5) Ridky, T.; Leis, J. J. Biol. Chem. 1995, 270, 29621. (6) Condra, J. H.; Schleif, W. A.; Blahy, O. M.; Gabryelski, L. J.; Graham, D. J.; Quintero, J. C.; Rhodes, A.; Robbins, H. L.; Roth, E.; Shivaprakash, M.; Titus, D.; Yang, T.; Teppler, H.; Squires, K. E.; Deutsch, P. J.; Emini, E. A. Nature 1995, 374, 569. (7) Ohtaka, H.; Muzammil, S.; Schon, A.; Velazquez-Campoy, A.; Vega, S.; Freire, E. Int. J. Biochem. Cell Biol. 2004, 36, 1787. (8) Lam, P. Y. S.; Jadhav, P. K.; Eyermann, C. J.; Hodge, C. N.; Ru, Y.; Bacheler, L. T.; Meek, J. L.; Otto, M. J.; Rayner, M. M.; Wong, Y. N.; Chang, C. H.; Weber, P. C.; Jackson, D. A.; Sharpe, T. R.; Ericksonviitanen, S. Science 1994, 263, 380. (9) Ala, P. J.; Huston, E. E.; Klabe, R. M.; McCabe, D. D.; Duke, J. L.; Rizzo, C. J.; Korant, B. D.; DeLoskey, R. J.; Lam, P. Y. S.; Hodge, C. N.; Chang, C. H. Biochemistry 1997, 36, 1573. (10) Krohn, A.; Redshaw, S.; Ritchie, J. C.; Graves, B. J.; Hatada, M. H. J. Med. Chem. 1991, 34, 3340. (11) Kempf, D. J.; Marsh, K. C.; Denissen, J. F.; McDonald, E.; Vasavanonda, S.; Flentge, C. A.; Green, B. E.; Fino, L.; Park, C. H.; Kong, X. P.; Wideburg, N. E.; Saldivar, A.; Ruiz, L.; Kati, W. M.; Sham, H. L.; Robins, T.; Stewart, K. D.; Hsu, A.; Plattner, J. J.; Leonard, J. M.; Norbeck, D. W. Proc. Natl. Acad. Sci. U.S.A. 1995, 92, 2484. (12) Ishima, R.; Freedberg, D. I.; Wang, Y. X.; Louis, J. M.; Torchia, D. A. Structure 1999, 7, 1047. (13) Freedberg, D. I.; Ishima, R.; Jacob, J.; Wang, Y. X.; Kustanovich, I.; Louis, J. M.; Torchia, D. A. Protein Sci. 2002, 11, 221. (14) Katoh, E.; Louis, J. M.; Yamazaki, T.; Gronenborn, A. M.; Torchia, D. A.; Ishima, R. Protein Sci. 2003, 12, 1376. (15) Scott, W. R. P.; Schiffer, C. A. Structure 2000, 8, 1259. (16) Meagher, K. L.; Carlson, H. A. Proteins: Struct., Funct., Bioinf. 2005, 58, 119. (17) Tóth, G.; Borics, A. J. Mol. Graphics Model. 2006, 24, 465. (18) Ishima, R.; Louis, J. M. Proteins: Struct., Funct., Bioinf. 2008, 70, 1408. (19) Collins, J. R.; Burt, S. K.; Erickson, J. W. Nat. Struct. Biol. 1995, 2, 334. (20) Hornak, V.; Okur, A.; Rizzo, R. C.; Simmerling, C. Proc. Natl. Acad. Sci. U.S.A. 2006, 103, 915. (21) Ding, F.; Layten, M.; Simmerling, C. J. Am. Chem. Soc. 2008, 130, 7184. (22) Hamelberg, D.; McCammon, J. A. J. Am. Chem. Soc. 2005, 127, 13778.
The Roles of Local Dynamics at the Flap Tips (23) Wiley, A. P.; Williams, S. L.; Essex, J. W. J. Chem. Theory Comput. 2009, 5, 1117. (24) Tozzini, V.; McCammon, J. A. Chem. Phys. Lett. 2005, 413, 123. (25) Tozzini, V.; Trylska, J.; Chang, C. E.; McCammon, J. A. J. Struct. Biol. 2007, 157, 606. (26) Chang, C. E.; Shen, T.; Trylska, J.; Tozzini, V.; McCammon, J. A. Biophys. J. 2006, 90, 3880. (27) Trylska, J.; Tozzini, V.; Chang, C. A.; McCammon, J. A. Biophys. J. 2007, 92, 4179. (28) Li, D. C.; Liu, M. S.; Ji, B. H.; Hwang, K.; Huang, Y. G. J. Chem. Phys. 2009, 130, 215102. (29) Li, D.; Ji, B.; Hwang, K.; Huang, Y. Int. J. Appl. Mech. 2009, 1, 113. (30) Berendsen, H. J. C.; Vanderspoel, D.; Vandrunen, R. Comput. Phys. Commun. 1995, 91, 43. (31) Lindahl, E.; Hess, B.; van der Spoel, D. J. Mol. Model. 2001, 7, 306. (32) Sorin, E. J.; Pande, V. S. Biophys. J. 2005, 88, 2472. (33) Junmei, W.; Wei, W.; Kollman, P. A.; Case, D. A. J. Mol. Graphics Model. 2006, 25, 247. (34) Wang, J. M.; Wolf, R. M.; Caldwell, J. W.; Kollman, P. A.; Case, D. A. J. Comput. Chem. 2004, 25, 1157. (35) Jakalian, A.; Bush, B. L.; Jack, D. B.; Bayly, C. I. J. Comput. Chem. 2000, 21, 132. (36) Case, D. A.; Cheatham, T. E.; Darden, T.; Gohlke, H.; Luo, R.; Merz, K. M.; Onufriev, A.; Simmerling, C.; Wang, B.; Woods, R. J. J. Comput. Chem. 2005, 26, 1668. (37) Jorgensen, W. L.; Chandrasekhar, J.; Madura, J. D.; Impey, R. W.; Klein, M. L. J. Chem. Phys. 1983, 79, 926. (38) Essmann, U.; Perera, L.; Berkowitz, M. L.; Darden, T.; Lee, H.; Pedersen, L. G. J. Chem. Phys. 1995, 103, 8577. (39) Berendsen, H. J. C.; Postma, J. P. M.; Vangunsteren, W. F.; Dinola, A.; Haak, J. R. J. Chem. Phys. 1984, 81, 3684. (40) Ryckaert, J. P.; Ciccotti, G.; Berendsen, H. J. C. J. Comput. Phys. 1977, 23, 327.
J. Phys. Chem. B, Vol. 114, No. 8, 2010 3069 (41) Humphrey, W.; Dalke, A.; Schulten, K. J. Mol. Graphics 1996, 14, 33. (42) Spinelli, S.; Liu, Q. Z.; Alzari, P. M.; Hirel, P. H.; Poljak, R. J. Biochimie 1991, 73, 1391. (43) Vriend, G. J. Mol. Graphics 1990, 8, 52. (44) Yamazaki, T.; Nicholson, L. K.; Torchia, D. A.; Wingfield, P.; Stahl, S. J.; Kaufman, J. D.; Eyermann, C. J.; Hodge, C. N.; Lam, P. Y. S.; Ru, Y.; Jadhav, P. K.; Chang, C. H.; Weber, P. C. J. Am. Chem. Soc. 1994, 116, 10791. (45) Ode, H.; Neya, S.; Hata, M.; Sugiura, W.; Hoshino, T. J. Am. Chem. Soc. 2006, 128, 7887. (46) Wittayanarakul, K.; Aruksakunwong, O.; Saen-oon, S.; Chantratita, W.; Parasuk, V.; Sompornpisut, P.; Hannongbua, S. Biophys. J. 2005, 88, 867. (47) Wittayanarakul, K.; Hannongbua, S.; Feig, M. J. Comput. Chem. 2008, 29, 673. (48) Perryman, A. L.; Lin, J. H.; McCammon, J. A. Protein Sci. 2004, 13, 1108. (49) Nicholson, L. K.; Yamazaki, T.; Torchia, D. A.; Grzesiek, S.; Bax, A.; Stahl, S. J.; Kaufman, J. D.; Wingfield, P. T.; Lam, P. Y. S.; Jadhav, P. K.; Hodge, C. N.; Domaille, P. J.; Chang, C. H. Nat. Struct. Biol. 1995, 2, 274. (50) Lipari, G.; Szabo, A. J. Am. Chem. Soc. 1982, 104, 4546. (51) Tjandra, N.; Wingfield, P.; Stahl, S.; Bax, A. J. Biomol. NMR 1996, 8, 273. (52) Duan, L. L.; Tong, Y.; Mei, Y.; Zhang, Q. G.; Zhang, J. Z. H. J. Chem. Phys. 2007, 127. (53) Markgren, P. O.; Schaal, W.; Hämäläinen, M.; Karlen, A.; Hallberg, A.; Samuelsson, B.; Danielson, U. H. J. Med. Chem. 2002, 45, 5430. (54) Martin, P.; Vickrey, J. F.; Proteasa, G.; Jimenez, Y. L.; Wawrzak, Z.; Winters, M. A.; Merigan, T. C.; Kovari, L. C. Structure 2005, 13, 1887. (55) Layten, M.; Hornak, V.; Simmerling, C. J. Am. Chem. Soc. 2006, 128, 13360. (56) Luzar, A.; Chandler, D. Nature 1996, 379, 55.
JP1005549
