Using Steered Molecular Dynamics to Predict and Assess Hsp70 Substrate-Binding Domain Mutants that Alter Prion Propagation

Abstract

Genetic screens using Saccharomyces cerevisiae have identified an array of cytosolic Hsp70 mutants that are impaired in the ability to propagate the yeast [PSI+] prion. The best characterized of these mutants is the Ssa1 L483W mutant (so-called SSA1-21), which is located in the substrate-binding domain of the protein. However, biochemical analysis of some of these Hsp70 mutants has so far failed to provide major insight into the specific functional changes in Hsp70 that cause prion impairment. In order to gain a better understanding of the mechanism of Hsp70 impairment of prions we have taken an in silico approach and focused on the Escherichia coli Hsp70 ortholog DnaK. Using steered molecular dynamics simulations (SMD) we demonstrate that DnaK variant L484W (analogous to SSA1-21) is predicted to bind substrate more avidly than wild-type DnaK due to an increase in numbers of hydrogen bonds and hydrophobic interactions between chaperone and peptide. Additionally the presence of the larger tryptophan side chain is predicted to cause a conformational change in the peptide-binding domain that physically impairs substrate dissociation. The DnaK L484W variant in combination with some SSA1-21 phenotypic second-site suppressor mutations exhibits chaperone-substrate interactions that are similar to wild-type protein and this provides a rationale for the phenotypic suppression that is observed. Our computational analysis fits well with previous yeast genetics studies regarding the functionality of the Ssa1-21 protein and provides further evidence suggesting that manipulation of the Hsp70 ATPase cycle to favor the ADP/substrate-bound form impairs prion propagation. Furthermore, we demonstrate how SMD can be used as a computational tool for predicting Hsp70 peptide-binding domain mutants that impair prion propagation.

Abstract

Genetic screens using Saccharomyces cerevisiae have identified an array of cytosolic Hsp70 mutants that are impaired in the ability to propagate the yeast [PSI+] prion. The best characterized of these mutants is the Ssa1 L483W mutant (so-called SSA1-21), which is located in the substrate-binding domain of the protein. However, biochemical analysis of some of these Hsp70 mutants has so far failed to provide major insight into the specific functional changes in Hsp70 that cause prion impairment. In order to gain a better understanding of the mechanism of Hsp70 impairment of prions we have taken an in silico approach and focused on the Escherichia coli Hsp70 ortholog DnaK. Using steered molecular dynamics simulations (SMD) we demonstrate that DnaK variant L484W (analogous to SSA1-21) is predicted to bind substrate more avidly than wild-type DnaK due to an increase in numbers of hydrogen bonds and hydrophobic interactions between chaperone and peptide. Additionally the presence of the larger tryptophan side chain is predicted to cause a conformational change in the peptide-binding domain that physically impairs substrate dissociation. The DnaK L484W variant in combination with some SSA1-21 phenotypic second-site suppressor mutations exhibits chaperone-substrate interactions that are similar to wild-type protein and this provides a rationale for the phenotypic suppression that is observed. Our computational analysis fits well with previous yeast genetics studies regarding the functionality of the Ssa1-21 protein and provides further evidence suggesting that manipulation of the Hsp70 ATPase cycle to favor the ADP/substrate-bound form impairs prion propagation. Furthermore, we demonstrate how SMD can be used as a computational tool for predicting Hsp70 peptide-binding domain mutants that impair prion propagation.

Author Summary

Direct non-covalent interactions between protein substrates and molecular chaperones play crucial roles in the protein folding process. [PSI+] is a prion of the yeast Saccharomyces cerevisiae, which is formed by mis-folding of the native Sup35 protein in a process analogous to formation of prions in mammals. While much genetic data exists showing a clear role for Hsp70 in prion propagation, biochemical data has yet to provide a clear link to how Hsp70 functions in prion propagation or how some Hsp70 mutants successfully impair in vivo propagation of prions. This paper employs a novel simulation method termed "steered molecular dynamics" to explore the different types and amounts of non-covalent interactions between wild type and mutated Hsp70, with a model substrate. Extrapolating the in silico data allowed us to decipher how a mutant Hsp70 impairs yeast prion propagation and allows us to predict other Hsp70 mutants that should behave in the same manner and to test these predictions in a yeast-based system. Our computational data shows that increasing the binding affinity of Hsp70 for substrate is one way of impairing prion propagation, a proposal that correlates very well with previous experimental genetic data.

Citation: Xu L, Hasin N, Shen M, He J, Xue Y, et al. (2013) Using Steered Molecular Dynamics to Predict and Assess Hsp70 Substrate-Binding Domain Mutants that Alter Prion Propagation. PLoS Comput Biol 9(1): e1002896. doi:10.1371/journal.pcbi.1002896

Editor: Giuseppe Legname, Scuola Internazionale Superiore di Studi Avanzati, Italy

Received: July 25, 2012; Accepted: December 8, 2012; Published: January 31, 2013

Copyright: © 2013 Xu et al. This is an open-access article distributed under the terms of the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original author and source are credited.

Funding: This work was supported by grants from the National Natural Science Foundation of China (No. 30970152) awarded to YS and a Science Foundation Ireland Research Frontiers award (RFP/07/BIC493) to GWJ. NH is a recipient of a NUI Maynooth John and Pat Hume postgraduate scholarship. SP was supported by the 973 Program (2012CB911000, 2013CB910700) and the National Natural Science Foundation of China (31110103914, 31070656). The funders had no role in study design, data collection and analysis, decision to publish, or preparation of the manuscript.

Competing interests: The authors have declared that no competing interests exist.

Introduction

Hsp70 is a highly conserved essential molecular chaperone that assists protein folding in general and protects cells from stress by preventing aggregation of stress-denatured proteins [1]. Hsp70 proteins consist of two functional domains: the N-terminal region comprises the nucleotide binding domain, and the C-terminal region the substrate-binding domain, which consists of a β-sheet rich substrate-binding cavity and a α-helical region that forms a lid over this cavity [2], [3]. N-terminal ATP hydrolysis induces a conformational change in the protein that closes the lid, trapping the substrate [3]. Additionally, binding of substrate causes a conformational change within the substrate-binding domain that stimulates ATP hydrolysis in proportion to its binding affinity, thus demonstrating a clear inter-domain communication mechanism for regulating Hsp70 function [4], [5].

In line with yeast prions being viewed as a protein-folding problem, altered abundance or function of a variety of molecular chaperones can affect their propagation [6]. Hsp104, a molecular chaperone that protects cells from exposure to stress by resolubilizing denatured protein from aggregates, is thought to facilitate yeast prion [PSI+] "seed" generation by breaking Sup35p aggregates into more numerous pieces [7]-[9]. Protein disaggregation by Hsp104 requires the assistance of Hsp70, which normally acts in [PSI+] "seed" propagation by dismantling polymers from Sup35 aggregates [10]. S. cerevisiae has four cytosolic Hsp70s of the Ssa subfamily, and expression of at least one of them is essential for growth [11]. Ssa1 is one of two constitutively expressed Ssa Hsp70s. Previous studies have identified and characterized a dominant mutant of Ssa1L483W (Ssa1-21) that impairs propagation of the yeast [PSI+] prion, weakening the phenotype caused by [PSI+] and causing a reduction in the number of [PSI+] "seeds" per cell [12], [13]. Despite the considerable impairment of [PSI+], Ssa1-21 has little effect on cell growth under optimal or stressful conditions, which therefore provides a convenient model to study yeast prion regulation in vivo [14], [15]. Previous in vivo studies suggested that Ssa1-21 has enhanced substrate binding compared with wild-type Ssa1 [14]. Due to the location of residue L483 in the substrate-binding domain, Ssa1-21 might bind more avidly to Sup35 aggregates and thus impair the ability Hsp70 to promote the disassembly of aggregates into smaller oligomers [14]. This hypothesis was further supported by observations using fluorescence recovery after photobleaching (FRAP) and semi-denaturing detergent-agarose gel electrophoresis (SDD-AGE) demonstrating that Sup35 aggregates are larger in SSA1-21 [PSI+] compared to wild type cells [13]. A recent in vitro study suggests that Ssa1-21 actually has weakened substrate-binding activity towards a small peptide substrate compared with wild-type Ssa1 [16]. The apparent discrepancy between the in vivo and in vitro activities of the Ssa1-21 protein may simply reflect the limitations of the biochemical assay used to mimic interaction with a prion substrate in vivo.

Molecular dynamics (MD) simulation is a classical method, which can provide detailed information on structural changes. However, applying such a technique to studying SBD- substrate interactions for Hsp70 would provide minimal information as the complete release of substrate from the SBD hydrophobic pocket would not be accomplished in the time scale of a conventional MD simulation. Steered molecular dynamics (SMD) is a way to imitate atomic force microscopy (AFM) experiments and can provide insight to the molecular mechanisms underlying the dissociation process of a protein and associated substrate [17], [18]. SMD simulations have been successfully used to investigate ligand-receptor interactions and the unfolding process of proteins [17]-[20]. We therefore applied the powerful SMD simulation technique to study the dissociation process of substrate from the SBD hydrophobic binding pocket of the highly conserved bacterial Hsp70 protein DnaK. DnaK is the most intensively studied Hsp70 protein with a high quality crystal structure available as well as detailed biochemical data [3]-[5], [21]. Corresponding structural information is severely lacking for Ssa1. Additionally, justified correlations have been made between DnaK biochemical activities (wild-type and mutant variants) and mutants obtained in the yeast ortholog Ssa1 [4], [22], [23]. Here we demonstrate an excellent correlation between previous experimental data for DnaK- substrate binding and that obtained through SMD simulations. We also use SMD simulations to characterize a DnaKL484W mutant that provides strong evidence suggesting that increased substrate binding is the reason behind the prion-impairing ability of the equivalent yeast Ssa1 mutant. Additionally we demonstrate that SMD can be used to identify and predict Hsp70 SBD mutants that will impair prion propagation in vivo. We further validate our approach of analyzing the well-characterized DnaK system by demonstrating that SMD simulations carried out on computer-modeled wild type and mutant Ssa1 proteins behave in a strikingly similar manner to their DnaK counterparts.

Results

MD equilibrium simulations for DnaKWT and DnaKL484W

The system temperature was coupled at 300K and the total energy fluctuated within less than 0.7% following short-time equilibration. The structural stability of models was assessed by calculation of the root-mean-square deviation (RMSD). To select equilibrated states as input structures for the SMD simulations, the overall RMSDs of Cα atoms of the protein were calculated as a function of the simulation time (Figure S1). The RMSD profiles indicate that the conformation of the proteins was equilibrated at approximately 2.5 ns and the values of RMSD for both DnaKWT and DnaKL484W were generally stabilized near 0.25 nm. The stable structures of each of the proteins provide appropriate starting points for subsequent SMD simulations and is also suggestive of similar structural stabilities for both proteins.

Determination of the pulling velocity and spring constant

In constant velocity SMD (cv-SMD), the value of the spring force varies significantly depending on what pulling velocity and spring constant were applied to the system [24]-[26]. Higher pulling velocities result in dis-equilibrium, which introduces significant error into the simulation results [25]. To mimic AFM experiments as much as possible, the pulling velocity should be as small as possible [27]. We carried out 12 different pulling velocity simulations (ranging from 0.001 nm/ps to 0.08 nm/ps) to identify a suitable velocity in this DnaK system. Figure S2A summarizes the results of the rupture forces caused by different pulling velocities. The rupture force was defined as the highest peak value of the pulling force. There was an apparent linear dependency between rupture force and pulling velocity when velocities were lower than 0.03 nm/ps. For subsequent SMD simulations a pulling velocity of 0.001 nm/ps was chosen, as the external harmonic potential has no obvious influence on the structure of the protein.

Another important consideration for SMD simulations is that the spring constant must be large enough to enable the local dissociation potential to be measured; however it must not be too large as this leads to difficulty in measurements due to high-levels of background noise [27]. In this study, we assessed six different spring constants (2000, 1000, 500, 300, 100, 50 kJ/mol/nm2). Figure S2B shows that under the same pulling velocity, the spring constants of 2000, 1000 and 500 kJ/mol/nm2 caused large fluctuations in rupture force, while spring constants of 100 and 50 kJ/mol/nm2 had minimal effect. A spring constant of 300 kJ/mol/nm2, however, produced a measurable and consistent level of rupture forces and was therefore chosen for subsequent SMD simulations.

Interaction energy analysis of the substrate dissociation process

The release of substrate from the cleft of a protein or enzyme is a dynamic process, which can be analyzed by measuring the interaction energy [27]. Following equilibration of the DnaK substrate bound systems and applying the appropriate previously derived parameters, we analyzed the interaction energy profiles for DnaKWT and DnaKL484W with substrate, during the simulated process of pulling substrate from the SBD cleft. Figure 1 shows a diagrammatic representation of the process and also shows the fluctuations in interaction energy for DnaKWT and DnaKL484W during the 2 ns SMD simulations. With the application of steering forces to the system dis-equilibrium is achieved by the first time-point. Following the time-course of the simulation, the interaction energy initially remained at a low level of about 1300 ps, which indicates the presence of strong interactions between peptide substrate and the DnaK SBD cleft during this dynamic procedure. However, as the substrate is eventually pulled away from the cleft, the energy curve quickly approaches zero. The energy curves for peptide dissociation can be roughly divided into two distinct phases: Phase I, when the peptide is located within the SBD binding site, which ranges from 0-1300 ps for DnaKWT and 0-1400 ps for DnaKL484W; and Phase II, when the peptide is released into the bulk solvent, which ranges from 1300-2000 ps for DnaKWT and 1400-2000 ps for DnaKL484W. The DnaKL484W-substrate complex clearly has a lower energy curve than the DnaKWT complex, especially in Phase I, which indicates a more stable DnaKL484W complex compared to DnaKWT. It therefore appears that although L484 is not located directly in the SBD cleft, this residue plays a key role in influencing DnaK substrate binding affinity.

[Figure omitted, see PDF]

Figure 1. The substrate release process of DnaK during 2 ns SMD.

Ribbon schematic representations of the DnaK SBD (393-507) with the peptide NRLLLTG, used in the SMD simulations. The simulated pulling of substrate from the binding cleft of DnaK was carried out over a 2000 ps time period. The DnaK SBD is yellow, substrate is red, and the residue 484 tryptophan is blue. The graph records changes in interaction energy during the course of the SMD simulation with DnaKL484W (red) and DnaKWT (black).

Hydrogen bond and hydrophobic interaction analyses during substrate release

Hydrogen bonds and hydrophobic interactions (or van der Waals) are important non-covalent interactions in stabilizing protein structures. Previous experimental research has indicated that the peptide interacts with the SBD of DnaK through van der Waals forces and seven intermolecular hydrogen bonds [3]; these interactions are summarized in Table 1. In our simulation studies, these two types of non-covalent interactions were the only ones detected in the system. Figure 2A shows the total number of hydrogen bonds formed between the peptide and DnaK as a function of time during the simulation. According to this curve, we can clearly determine that during the simulation the total number of intermolecular hydrogen bond between DnaKL484W and substrate is greater than that of DnaKWT. Moreover, in our SMDs five out of the seven major hydrogen bonds identified by Zhu et al. [3] contribute significantly to stabilizing the DnaKWT complex while four contribute to DnaKL484W complex stabilization (Table 1). However, while having a reduced number of hydrogen bond interactions compared to wild type, those formed between DnaKL484W and substrate are in existence for a much longer time period during the time course of the simulation (Table 1). Three of these, Ala429(N)-Leu904(O), Ser427(N)-Arg902(O) and Ser427(O)-Leu904(N), are the most significant intermolecular hydrogen bonds influencing the dissociation processes. In addition, it was found that the lifetimes of one extra hydrogen bond (Thr428(OG1)-Thr906(OG1)) in the DnaKWT complex, and three extra hydrogen bonds (Thr409(OG1)-Asn901(N), Thr428(OG1)-Thr906(OG1), Ala435(O)-Gly907(N)) in the DnaKL484W complex, are longer by 50% in our simulations. This suggests these bonds may also play important roles in influencing the dissociation process (Table 1).

[Figure omitted, see PDF]

Figure 2. Changes in interaction forces as a function of simulation time.

(A) Temporal evolution of the number of hydrogen bonds number formed between the SBD of DnaK and peptide under the application of force for DnaKL484W (red) and DnaKWT (black) complex; (B) Variations in the numbers of direct hydrophobic interactions of DnaKL484W (red) and DnaKWT (black) complexes in the SMD simulation, respectively. The calculation of interactions is recorded every 100 ps.

Table 1. Comparison of peptide-substrate interactions identified in SMD simulations with previous experimental data.

In this study the number of van der Waals contacts of nonpolar atoms represents the hydrophobic interaction strength [28]. As shown in Figure 2B, the average number of hydrophobic interactions of the DnaKL484W complex is obviously greater than that of DnaKWT. Moreover, the rupture behavior of the hydrophobic interactions at the binding site in DnaKL484W complex occurred later than that in the DnaKWT complex, which is similar for the rupture behavior of the hydrogen bond network (Figures 2A and 2B). Moreover, in our SMDs most of hydrophobic interactions identified by Zhu et al. [3] contribute significantly to stabilizing the DnaKWT and DnaKL484W complexes (Table 1).

To analyze in more detail the interactions between the substrate and the hydrophobic pocket, specific structural snapshots were taken at the largest rupture force point (DnaKWT at 1300 ps, DnaKL484W at 1400 ps). Hydrogen bonds and hydrophobic interactions between the substrate and binding pocket were analyzed using the LIGPLOT program. As shown in Figure 3 the three most important intermolecular hydrogen bonds, Ala429(N)- Leu904(O), Ser427(N)-Arg902(O) and Ser427(O)-Leu904(N) , were still formed in DnaKL484W at this moment, whereas only two of them, Ala429(N)-Leu904(O) and Ser427(N)-Arg902(O), were formed in DnaKWT. Additionally, at this point there are six hydrogen bonds in the DnaKL484W complex and only four in the case of DnaKWT. Furthermore, 101 hydrophobic interactions were detected between residues of the SBD (Thr403, Met404, Val425, Phe426, Ser427, Thr428, Ala429, Val436, Gln471) and substrate in DnaKL484W, especially Thr403, Met404, and Ser427, compared to only 81 hydrophobic interactions in DnaKWT (formed between Met404, Gln424, Phe426, Ser427, Thr428, Ala429, Glu430, Val436 and substrate).

[Figure omitted, see PDF]

Figure 3. Interactions at the largest rupture force.

Schematic diagrams showing hydrogen bonds and hydrophobic interactions between substrate and substrate-binding domain complexes for DnaKWT (A) and DnaKL484W (B) at the largest rupture force. The hydrophobic interactions (with distances of less than 3.5 Å) are shown as spiked spheres and dashed lines indicate hydrogen bonds.

The detailed interaction analysis from this in silico work, as well as interaction analysis from previous experimental data is summarized in Table 1. The most significant interactions between SBD residues and substrate identified in this study are highly consistent with previous experimental data [3]. This fact validates the accuracy of the simulation and robustness of the data obtained.

Conformational changes in the SBD in the L484W mutant

L484W is a mutation within the C-terminal region of DnaK. Previous work suggests that the equivalent mutation in yeast cytosolic Hsp70 causes increased substrate binding [14], a proposal that is clearly reflected in our binding analysis. Although L484W substitution is within the SBD, it is distant from the substrate-binding pocket. So its effect on the substrate-binding affinity may be indirect perhaps through conformational changes induced within the SBD. Previous studies have shown that the arch structure formed by residues Met404 and Ala429 of the two substrate contacting loops in the SBD of DnaK contributes significantly to peptide binding [3], [5], [21]. These two residues were also highlighted in both DnaKWT and DnaKL484W in our SMD simulations (Table 1, Figure 3). Using the VMD program we assessed potential structural changes within the SBD between DnaKWT and DnaKL484W. The final structure (following MD simulations) of DnaKL484W was superimposed with the corresponding final structure of DnaKWT for the Cα (Figure 4A). We found that DnaKL484W can more readily extend Ala429 to Met404 and then form a tighter arch than DnaKWT (Figures 4B and 4C). It seems therefore that L484W, a tryptophan mutant in β7, somehow induces long-range conformational changes that result in the formation of a tighter arch structure that may impede efficient substrate release (Figure 4A).

[Figure omitted, see PDF]

Figure 4. Superimposition of DnaKWT and DnaKL483W structures.

(A) Superimposition of the Cα atoms of DnaKL484W (red) and DnaKWT (black) structures at the end of the MD simulation. (B) and (C) show superposition of the arch structure snapshots for DnaKWT and DnaKL484W, respectively.

Analysis of other Leucine 484 mutations

In order to explore further the possible regulation mechanism of the L484W variant, we performed more MD and SMD simulations using different amino acid substitutions at L484, namely DnaKL484A and DnaKL484H. Alanine has a similar hydrophobicity as tryptophan but has a much smaller side chain while histidine has a similar size of side chain as tryptophan but a much lower hydrophobicity. With these two "control" mutants, we can investigate whether the hydrophobicity and/or side chain structure may play important roles in the L484W variant. Compared with DnaKWT, DnaKL484W and DnaKL484H, the DnaKL484A mutant shows significant fluctuations and a higher RMSD value during the simulations (the RMSD averages of Cα atoms are shown in Table 2). This suggests that the DnaKL484A variant may be more structurally unstable than DnaKL484H. In addition to differences in distance of arch structure residues, we also detected major differences in interaction energy, numbers of hydrogen bonds and hydrophobic interactions between the DnaKL484A and DnaKL484H variants (Table 2). In summary, DnaKL484A showed substrate interaction characteristics similar to DnaKWT, while DnaKL484H exhibited characteristics akin to DnaKL484W (Table 2). These findings suggest that the structural and functional changes caused by L484W are likely due to the influence of the side chain rather than changes in hydrophobicity.

Table 2. Effects on substrate-binding interaction of different DnaK mutations.

Increased SBD-substrate interactions are responsible for L484W effects

The results from the in silico structural analysis of DnaKL484A and DnaKL484H affords the opportunity to formally assess the functional impact of increased substrate binding upon DnaK. The interaction profiles would predict that DnaKL484A would behave in a similar manner to DnaKWT while DnaKL484H would behave similarly to DnaKL484W (Table 2).

We therefore exploited the yeast [PSI+] system to test this hypothesis. Previous studies identified the Ssa1L483W mutant (equivalent to DnaKL484W) as impairing propagation of the yeast prion [PSI+] [12] and genetics studies suggest that this mutant binds more avidly to substrate [14]. To assess whether the SMD studies performed on DnaK could potentially be translated into meaningful comparisons with the yeast Ssa1 protein we constructed a molecular model of Ssa1 in complex with the same substrate as for DnaK, and performed the same SMD analysis on wild-type and L483 mutants as was carried out for DnaK (Text S1). Strikingly, the SMD analysis for Ssa1 is in excellent agreement with the SMD analysis for DnaK. Both Ssa1WT and Ssa1L483W display similar stability profiles (Figure S3) and Ssa1L483W displayed an overall tighter affinity for substrate (Figure S4). Additionally Ssa1L483W displayed an overall higher number of both hydrogen bonds and hydrophobic interactions over the time course of the SMD and the interaction profiles of the corresponding L483A and L483H mutants correlated well with their DnaK counterparts (Figures S5, Table S1).

We therefore proceeded with in vivo studies in yeast and constructed Ssa1L483A and Ssa1L483H variants and tested their ability to propagate [PSI+] when they were the sole cytosolic Hsp70 in the cell. Figure 5 demonstrates that similar to Ssa1WT the Ssa1L483A mutant can efficiently propagate [PSI+] while similar to Ssa1L483W the Ssa1L483H variant has impaired [PSI+] propagation; precisely what was predicted from the in silico DnaK studies. Additionally, we found that unlike Ssa1L483W the Ssa1L483H mutant could propagate [PSI+], albeit at a much-reduced level compared to Ssa1WT (Figure 5, -ADE 4 day panel), which may be a reflection of the wild-type level of hydrogen bonds between inner loop residues (Table S1) preventing the full impairment of the prion exhibited by L483W. The experimental confirmation of the in silico predictions demonstrates that MD and SMD simulations can be used to predict and design Hsp70 SBD mutants that may impair prion propagation in vivo.

[Figure omitted, see PDF]

Figure 5. Assessment of [PSI+] prion propagation in L483 mutants.

Growth on adenine and color on YPD is indicative of yeast prion status. [PSI+] and [psi-] controls are in the presence of wild-type Ssa1. Strains are only expressing one type of cytosolic Hsp70, indicated by mutant type.

Finally, previous experimental studies had identified second-site suppressor mutations of Ssa1L483W that could restore [PSI+] prion propagation [14]. Mutating conserved L483W second-site suppressors in DnaK provides the ideal opportunity for in silico testing for effects on SBD-substrate interactions when present in combination with DnaKL484W. We therefore assessed the structural consequences of introducing either of two mutations, G406S or R445K, in combination with L484W. MD and SMD simulations demonstrate that when either of these mutations is placed in combination with L484W there is a restoration of interaction energy and a reduction in the numbers of hydrogen bonds and hydrophobic interactions to levels similar to wild type (Table 2). This finding further implies that the major cause of functional and phenotypic changes in the Hsp70 L483/4W mutation is due to increased SBD-substrate interaction.

Discussion

In this study we have taken an in silico approach to studying interactions between substrate and the SBD of the evolutionary conserved bacterial Hsp70 DnaK. Applying SMD simulations to DnaKWT and the DnaKL484W mutant has provided detailed structural information suggesting the L484W mutant has enhanced substrate-binding activity. Our structural results demonstrate that (i) the SBDs of DnaKL484W and DnaKWT have similar structural stability but the interaction energy of DnaKL484W with substrate is stronger, (ii) the intermolecular hydrogen bonds and hydrophobic interactions between substrate and DnaKL484W during the dissociation process is increased compared to DnaKWT, (iii) DnaKL484W can more efficiently form the critical Ala429 to Met404 interaction that results in formation of a tighter arch to prevent substrate release compared with DnaKWT, (iv) the side chain structure, not hydrophobicity, play an important role in DnaKL484W-induced functional changes. The simulation results obtained in this study compare well to previous experimental studies addressing DnaK-substrate interactions [3] and also provide further evidence to suggest that the yeast Ssa1L483W mutant and the equivalent DnaK mutant do indeed bind more avidly to substrate, as previously proposed [14].

Further validation of the use of SMD to study DnaK-substrate interactions comes from the analysis of the DnaKL484A and DnaKL484H mutants. Of particular note was the interaction profile of DnaKL484H, which clearly indicates that this mutant should have similar functional properties as the DnaKL484W mutant. A key goal of this work was to provide insight into the functional changes in the yeast SSA1-21 mutant that allow prion impairment. Through modeling wild-type Ssa1 and SBD mutants, and carrying out SMD simulations we clearly demonstrate that for the L483/4 residue of Ssa1/DnaK, the predicted functional changes that occur through mutation of this site are mimicked between the yeast and bacterial Hsp70 orthologs. Making the transition into the yeast system and demonstrating that when present in Ssa1 these two mutants significantly impair [PSI+] prion propagation, is a major finding that suggests SMD can be used to design Hsp70 SBD mutants that would be predicted to impair prion propagation. Such an application may be important for the possible future design of small therapeutic molecules that may target the SBD of Hsp70. Given the importance of Hsp70 in the protein folding process and the array of human protein-folding diseases, SMD could aid significantly in the early stages of future drug design and development for these diseases.

Previous genetic studies have shown that the Ssa1-21 mutant causes a reduction in seed/propagon numbers in vivo [12], [13]. The most likely cause for prion impairment was proposed to be through an increased affinity for prion substrate exhibited by Ssa1-21, which would interfere with the prion replication cycle [13]-[15]. Strong genetic support exists for this proposal in that an array of second-site mutations have been isolated that can suppress the ability of Ssa1-21 to impair prion propagation and biochemical assessment of corresponding mutations in DnaK has previously shown that such mutations reduce substrate-binding affinity. Therefore it seems that mutations predicted to reduce substrate-binding efficiency suppress the prion impairing ability of Ssa1-21, thus suggesting that Ssa1-21 impairment of [PSI+] is dependent on an increased substrate-binding activity of this mutant [14]. Additionally, deletion of co-chaperones known to promote substrate binding can also suppress the prion impairing ability of Ssa1-21 [15]. Further support for this hypothesis comes from our in silico studies on DnaK which not only predict increased substrate-binding affinity for the L484W variant but also predict a return to almost wild-type interaction profiles when L484W is in combination with either G406S or R445K. Both of these mutations were shown to suppress prion impairment exhibited by Ssa1L483W [14].

A recent biochemical study found that compared to wild type the Ssa1-21 protein possesses hyperactive intrinsic ATPase activity and a reduced capacity to bind a peptide substrate [16]. This data is somewhat at odds with the in vivo activity of Ssa1-21 since the mutant protein is able to maintain cell growth similar to wild type and as yet no stress-related phenotypes have been associated with this mutant. Similarly, our in silico data is in excellent agreement with in vivo data and hypothesis of Ssa1-21 prion impairment [13]-[15] but somewhat at odds with the biochemical data regarding substrate-binding efficiency [16]. There are a number of complicating factors that could account for this apparent discrepancy. Firstly, the peptide substrates used in the two studies are different. Biochemical peptide-binding assays utilized the A7 peptide (Arg-Arg-Leu-Glu-Asp-Ala-Glu-Tyr-Ala-Ala-Arg-Gly)whereas our in silico study relied on the substrate used in the crystallization of DnaK SBD (Asn-Arg-Leu-Leu-Leu-Thr-Gly). These peptide substrates are very different in amino acid composition and in length and therefore complicate any direct comparison between the in silico and in vitro data. Secondly, the biochemical analysis was carried out using full-length Ssa1 protein in the ADP-bound state, while the in silico analysis used the isolated SBD of DnaK for which high-quality structural data is available. Hence, it is difficult to draw direct comparison between current biochemical data for Ssa1-21 and our computational data. Our data clearly predicts increased binding affinities for both DnaKL484W and Ssa1L483W, which is supported by previous genetic studies and our own new in vivo data. Future biochemical analysis of these proteins should utilize a variety of peptides when assessing substrate-binding affinities and should also assess full-length and truncated versions of the chaperones, in order to understand the relationship between regulation of Ssa1 function in vivo and its biochemical properties in vitro.

In conclusion, we have demonstrated that SMD simulations can be used to study the process of substrate dissociation for DnaK and that in some cases this may apply to other Hsp70 orthologous. Of particular note is that SMD studies of Hsp70 may be used to predict and design new mutants that impair in vivo propagation of prions.

Materials and Methods

Model construction

In this study we used the solution structure of the substrate-binding domain from bacterial Hsp70 chaperone protein DnaK (393-507) in complex with the peptide NRLLLTG, obtained from the RCSD protein data bank (PDB code: 1Q5L) [21] as representing wild type protein. Appropriate models representing mutants of the L484 residue and subsequent second-site suppressor mutations, were constructed based upon the wild-type DnaK SBD structure using Swiss-Pdb Viewer [29] and subjected to Procheck analysis.

MD and SMD simulations

Simulations were carried out using the GROMACS 4.0.7 software package [30] with constant number, pressure, temperature and periodic boundary conditions. The GROMOS96 43a1 [31] force field was applied in all simulations. The models were immersed in the cubic boxes filled with water molecules with a distance between peptides and box edges of at least 10 Å. Each cubic box was filled with approximately 10315 simple point-charge water molecules. The linear constraint solver (LINCS) method [32] was used to constrain bond lengths, allowing an integration step of 2 fs. Electrostatic interactions were calculated using the particle mesh Ewald algorithm with a distance cut off of 9 Å [33], [34]. Lennard-Jones potential was used for describing the short-range attractive and repulsive dispersion interactions with a 10 Å cut off. The temperature and pressure were coupled using V-rescale and Parrinello-Rahman algorithm respectively. The solvated and neutralized systems were then energy-minimized for 2,000 steps using the Steepest Descent algorithm to remove bad contacts. The backbone atoms of the structure were then fixed, while the side chains and solvent were allowed to move unrestrained for 80 ps. After equilibration several independent simulations were carried out over a 5 ns period at 300 K, pH7.0 and 1bar pressure. The coordinate trajectories were saved every 1 ps for subsequent data analysis.

Following 5 ns simulations models in the equilibrated state were used as input structures for the COM pulling method of GROMACS, to visualize pulling of the substrate through the substrate-binding domain cleft (Figure 1). In this process, a soft elastic spring of 300 kJ mol-1 nm-2 was associated with the substrate and applied in the direction along the cleft with a constant velocity. The substrate-binding domain was considered as the static reference in the pulling simulation. The trajectories were saved every 0.1 ps, and steering forces were recorded every 10 fs. Multiple SMD experiments were employed to each system for 2 ns.

Analysis of trajectories

Analyses were performed using applications within the GROMACS package. All images were produced using the VMD [35]. The RMSD were calculated for all alpha-carbon atoms with reference to the first frame of the trajectories. Hydrogen bonds were calculated by a distance cutoff of 3.5 Å and an angle cutoff of 300°. Hydrophobic interactions were analyzed using Ligplot [36].

Yeast strains, media and plasmids

The yeast strain used in this study was G400-1C (MATa ade2-1 SUQ5 kar1-1 his3 leu2 lys2 trp1 ura3 ssa1::KanMX, ssa2::HIS3, ssa3::TRP1, ssa4::URA3-1f/pRDW10 [PSI+]) [14]. The plasmid used was pC210, which is a LEU2 based plasmid that contains the SSA1 coding sequence under control of the constitutively active SSA2 promoter [37]. Standard yeast media was used as previously described [38].

Site directed mutagenesis and yeast genetics

Mutants of L483 (either W/H/A) were introduced into plasmid pC210 by site-directed mutagenesis (SDM) using the Quickchange kit (Stratagene) as recommended by the manufacturer. For L48W the primers L483W-F, GACTCTAACGGTATTTGGAATGTTTCCGCC and L483W-R GGCGGAAACATTCCAAATACCGTTAGAGTC were used, for L483H the primers L483H-F, GACTCTAACGGTATTCACAATGTTTCCGCC and L483H-R, GGCGGAAACATTGTGAATACCGTTAGAGTC, and for introducing L483A the primers L483A-F, GACTCTAACGGTATTGCGAATGTTTCCGCC and L483A-R, GGCGGAAACATTCGCAATACCGTTAGAGTC were used. Underlined codons indicate the SDM change. Ssa1 L483 mutants were assessed for [PSI+] prion-impairing effects using the plasmid shuffle technique. Strain G400-1C was transformed with pC210 containing appropriate SSA1 L483 mutation. Transformed cells were selected on medium lacking leucine. Transformants were re-streaked onto fresh medium and then replica plated onto medium containing limiting amounts of adenine and also 5-fluro-orotic acid (5-FOA), a chemical that selects against URA+ cells and hence against the presence of the pRDW10 plasmid. Strains were then assessed for effects on the prion by using the well-defined red/white color assay based on the provision of adenine [39]. Briefly, the presence of [PSI+] (the non-functional aggregated form of Sup35p) and SUQ5 in strain G400-1c causes efficient translation read through of the ochre mutation in the ade2-1 allele. Non-suppressed ade2-1 mutants are Ade- and are red when grown on medium containing limiting amounts of adenine due to the accumulation of a pigmented substrate of Ade2p. Partial suppression of ade2-1 by [PSI+] allows growth without adenine and eliminates the pigmentation [40]. The presence or absence of [PSI+] was confirmed by mating experiments [39].

Supporting Information

Figure S1.

The Cα RMSD of the DnaKL484W (red) and DnaKWT (black) as a function of simulation time.

(JPG)

Figure S2.

Assessment of pulling velocity and spring constant in DnaK complex. (A) Computed rupture forces as a function of pulling velocity Vpull. The error bars give an estimated uncertainty. The dashed line shows a linear fit to the computed forces for Vpull less than 0.03 nm/ps. (B) Influence of different spring constants on the steering force of the DnaK complex.

(JPG)

Figure S3.

The Cα RMSD of the Ssa1L483W (red) and Ssa1WT (black) as a function of simulation time.

(JPG)

Figure S4.

Time dependence of the interaction energy between peptide and Ssa1 cleft for Ssa1L483W (red) and Ssa1WT (black).

(JPG)

Figure S5.

Changes in interaction forces as a function of simulation time. (A), Temporal evolution of the hydrogen bonds number formed between the SBD of Ssa1 and peptide under the application of force for Ssa1L483W (red) and Ssa1WT (black) complexes. (B) Variation in the numbers of direct hydrophobic interactions of Ssa1L483W (red) and Ssa1WT (black) complexes in the SMD simulations. The calculation of interactions is recorded every 100 ps.

(JPG)

Table S1.

Effects on substrate-binding interaction during the SMD simulation for different Ssa1 and different L483 mutant derivatives.

(PDF)

Text S1.

Construction of model for substrate-bound Ssa1.

(PDF)

Acknowledgments

We are grateful to Yu Wang for technical assistance in MD simulations and to Drs Hong Zhang and Weibin Gong for comments on the manuscript.

Author Contributions

Conceived and designed the experiments: LX YS GWJ. Performed the experiments: LX NH MS JH XZ. Analyzed the data: LX YX SP YS GWJ. Contributed reagents/materials/analysis tools: YS GWJ. Wrote the paper: LX SP YS GWJ.

References

Hartl FU, Hayer-Hartl M (2002) Molecular chaperones in the cytosol: from nascent chain to folded protein. Science 295: 1852-1858. doi: 10.1126/science.1068408. Find this article online

Morshauser RC, Wang H, Flynn GC, Zuiderweg ERP (1995) The peptide-binding domain of the chaperone protein HSC70 has an unusual secondary structure topology. Biochemistry 34: 6261-6266. doi: 10.1021/bi00019a001. Find this article online

Zhu X, Zhao X, Burkholder WF, Gragerov A, Ogata CM, et al. (1996) Structural analysis of substrate binding by the molecular chaperone DnaK. Science 272: 1606-1614. doi: 10.1126/science.272.5268.1606. Find this article online

Mayer MP, Schroder H, Rüdiger S, Paal K, Laufen T, et al. (2000) Multistep mechanism of substrate binding determines chaperone activity of Hsp70. Nat Struct Biol 7: 586-593. Find this article online

Pellecchia M, Montgomery DL, Stevens SY, Vander Kooi CW, Feng H, et al. (2000) Structural insights into substrate binding by the molecular chaperone DnaK. Nat Struct Biol 7: 298-303. Find this article online

Jones GW, Tuite MF (2005) Chaperoning prions: the cellular machinery for propagating an infectious protein? BioEssays 27: 823-832. doi: 10.1002/bies.20267. Find this article online

Chernoff YO, Lindquist SL, Ono B, Inge-Vechtomov SG, Liebman SW (1995) Role of the chaperone protein Hsp104 in propagation of the yeast prion-like factor [PSI+]. Science 268: 880-884. Find this article online

Kryndushkin DS, Alexandrow IM, Ter-Avanesyan MD, Kushnirov VV (2003) Yeast [PSI+] prion aggregates are formed by small Sup35 polymers fragmented by Hsp104. J Biol Chem 278: 49636-49643. Find this article online

Shorter J, Lindquist S (2004) Hsp104 catalyzes formation and elimination of self-replicating Sup35 Prion conformers. Science 304: 1793-1794. doi: 10.1126/science.1098007. Find this article online

Glover JR, Lindquist S (1998) Hsp104, Hsp70, and Hsp40: a novel chaperone system that rescues previously aggregated proteins. Cell 94: 73-82. Find this article online

Werner-Washburne M, Stone DE, Craig EA (1987) Complex interactions among members of an essential subfamily of Hsp70 genes in Saccharomyces cerevisiae. Mol Cell Biol 7: 2568-2577. Find this article online

Jung G, Jones G, Wegrzyn RD, Masison DC (2000) A role for cytosolic Hsp70 in yeast [PSI+] prion propagation and [PSI+] as a cellular stress. Genetics 156: 559-570. Find this article online

Song Y, Wu Y, Jung G, Tutar Y, Eisenberg E, et al. (2005) Role for Hsp70 chaperone in Saccharomyces cerevisiae prion seed replication. Eukaryotic Cell 4: 289-297. Find this article online

Jones GW, Masison DC (2003) Saccharomyces cerevisiae Hsp70 mutations affect [PSI+] prion propagation and cell growth differently and implicate Hsp40 and tetratricopeptide repeat cochaperones in impairment of [PSI+]. Genetics 163: 495-506. Find this article online

Jones G, Song Y, Chung S, Masison DC (2004) Propagation of Saccharomyces cerevisiae [PSI+] prion is impaired factors that regulate Hsp70 substrate binding. Mol Cel Bio 24: 3928-3937. Find this article online

Needham PG, Masison DC (2008) Prion-impairing mutations in Hsp70 chaperone Ssa1: effects on ATPase and chaperone activities. Arch Biochem Biophys 478(2): 167-174. doi: 10.1016/j.abb.2008.07.023. Find this article online

Kalikka J, Akola J (2011) Steered molecular dynamics simulations of ligand-receptor interaction in lipocalins. Eur Biophys J 40: 181-194. doi: 10.1007/s00249-010-0638-3. Find this article online

González A, Perez-Acle T, Leonardo Pardo, Deupi X (2011) Molecular basis of ligand dissociation in β-Adrenergic receptors. PLoS ONE 6(9): e23815 doi:10.1371/journal.pone.002381.5.

Hytönen VP, Vogel V (2008) How force might activate talin's vinculin binding sites: SMD reveals a structural mechanism. PLoS Comput Biol 4(2): e24 doi:10.1371/journal.pcbi.0040024.

Shen M, Guan J, Xu L, Yu Y, He J, et al. (2012) Steered molecular dynamics simulations on the binding of the appendant structure and helix-β2 in domain-swapped human cystatin C dimer. J Biomol Struct Dyn 30(6): 652-661. doi: 10.1080/07391102.2012.689698. Find this article online

Stevens SY, Cai S, Pellecchia M, Zuiderweg ERP (2003) The solution structure of the bacterial Hsp70 chaperone protein domain DnaK (393-507) in complex with the peptide NRLLLTG. Protein Science 12: 2588-2596. doi: 10.1110/ps.03269103. Find this article online

Montgomery DL, Morimoto RI, Gierasch LM (1999) Mutations in the substrate binding domain of the Escherichia coli 70 kDa molecular chaperone, DnaK, which alter substrate affinity or interdomain coupling. J Mol Biol 286: 915-932. Find this article online

Pfund C, Huang P, Lopez-Hoyo N, Craig EA (2001) Divergent functional properties of the ribosome-associated molecular chaperone Ssb compared with other Hsp70s. Mol Biol Cell 12: 3773-3782. Find this article online

Izrailev S, Stepaniants S, Balsera M, Oono Y, Schulten K (1997) Molecular dynamics study of unbinding of the avidin-biotin complex. Biophys J 72(4): 1568-1581. doi: 10.1016/S0006-3495(97)78804-0. Find this article online

Isralewitz B, Gao M, Schulten K (2001) Steered molecular dynamics and mechanical functions of proteins. Curr Opin Struct Biol 11: 224-230. doi: 10.1016/S0959-440X(00)00194-9. Find this article online

Santos SM, Costa PJ, Lankshear MD, Beer PD, Felix V (2010) Molecuar dynamics study of a Heteroditopic-Calix[4]diquinone-assisted transfer of KCl and dopamine through a water-chloroform Liquid-Liquid interface. J Phys Chem 114(34): 11173-11180. doi: 10.1021/jp100724e. Find this article online

Zhang J, Zheng Q, Zhang H (2010) Unbinding of glucose from human pulmonary surfactant protein D studied by steered molecular dynamics simulations. Chemical Physics Letters 484: 338-343. doi: 10.1016/j.cplett.2009.12.022. Find this article online

Shen L, Shen J, Luo X, Cheng F, Xu Y, et al. (2003) Steered molecular dynamics simulation on the binding of NNRTI to HIV-1 RT. J Biophysical 84: 3547-3563. doi: 10.1016/S0006-3495(03)75088-7. Find this article online

Guex N, Peitsch MC (1997) SWISS-MODEL and the Swiss-Pdb Viewer: an environment for comparative protein modeling. Electrophoresis 18: 2714-2723. doi: 10.1002/elps.1150181505. Find this article online

Van Der Spoel D, Lindahl E, Hess B, Groenhof G, Mark AE, et al. (2005) GROMACS: fast, flexible, and free. J Comput Chem 26: 1701-1718. doi: 10.1002/jcc.20291. Find this article online

Scott WRP, Hünenberger PH, Tironi IG, Mark AE, Billeter SR, et al. (1999) The GROMOS biomolecular simulation program package. J Phys Chem A 103: 3596-3607. doi: 10.1021/jp984217f. Find this article online

Hess B, Bekker H, Berendsen HJC, Fraaije J (1997) LINCS: a linear constraint solver for molecular simulations. J Comput Chem 18: 1463-1472. Find this article online

Darden T, York D, Pedersen L (1993) Particle mesh Ewald: an N-log(N) method for Ewald sums in large systems. J Chem Phys 98: 10089-10092. doi: 10.1063/1.464397. Find this article online

Essmann U, Perera L, Berkowitz ML, Darden T, Lee H, et al. (1995) A smooth particle mesh Ewald potential. J Chem Phys 103: 8577-8592. Find this article online

Humphrey W, Dalke A, Schulten K (1996) VMD: visual molecular dynamics. J Mol Graphics 14: 33-38. doi: 10.1016/0263-7855(96)00018-5. Find this article online

Wallace AC, Laskowski RA, Thornton JM (1995) LIGPLOT: a program to generate schematic diagrams of protein-ligand interactions. Protein Engineering 8: 127-134. doi: 10.1093/protein/8.2.127. Find this article online

Schwimmer C, Masison DC (2002) Antagonistic interactions between yeast [PSI+] and [URE3] prions and curing of [URE3] by Hsp70 protein chaperone Ssa1p but not by Ssa2p. Molecular and cellular biology 22: 3590-3598. Find this article online

Jones GW, Song Y, Masison DC (2003) Deletion of the Hsp70 chaperone gene SSB causes hypersensitivity to guanidine toxicity and curing of the [PSI+] prion by increasing guanidine uptake in yeast. Mol Gen Genomics 269: 304-311. Find this article online

Loovers HM, Guinan E, Jones GW (2007) Importance of the Hsp70 ATPase domain in yeast prion propagation. Genetics 175: 621-630. doi: 10.1534/genetics.106.066019. Find this article online

Cox BS (1965) [PSI+], a cytoplasmic suppressor of super-suppressors in yeast. Heredity 20: 505-521. Find this article online