Skip to main content
Advertisement
  • Loading metrics

Two-step mechanism of J-domain action in driving Hsp70 function

  • Bartlomiej Tomiczek ,

    Contributed equally to this work with: Bartlomiej Tomiczek, Wojciech Delewski, Lukasz Nierzwicki

    Roles Conceptualization, Data curation, Formal analysis, Investigation, Methodology, Resources, Software, Validation, Visualization, Writing – original draft, Writing – review & editing

    Affiliation Intercollegiate Faculty of Biotechnology, University of Gdansk and Medical University of Gdansk, Gdansk, Poland

  • Wojciech Delewski ,

    Contributed equally to this work with: Bartlomiej Tomiczek, Wojciech Delewski, Lukasz Nierzwicki

    Roles Investigation, Methodology, Writing – review & editing

    Affiliations Intercollegiate Faculty of Biotechnology, University of Gdansk and Medical University of Gdansk, Gdansk, Poland, Department of Biochemistry, University of Wisconsin-Madison, Madison, Wisconsin, United States of America

  • Lukasz Nierzwicki ,

    Contributed equally to this work with: Bartlomiej Tomiczek, Wojciech Delewski, Lukasz Nierzwicki

    Roles Conceptualization, Data curation, Formal analysis, Investigation, Methodology, Software, Validation, Visualization, Writing – original draft, Writing – review & editing

    Affiliation Department of Physical Chemistry, Gdansk University of Technology, Gdansk, Poland

  • Milena Stolarska,

    Roles Data curation, Investigation, Methodology, Software, Visualization, Writing – review & editing

    Affiliation Intercollegiate Faculty of Biotechnology, University of Gdansk and Medical University of Gdansk, Gdansk, Poland

  • Igor Grochowina,

    Roles Data curation, Investigation, Methodology, Resources, Visualization, Writing – review & editing

    Affiliation Intercollegiate Faculty of Biotechnology, University of Gdansk and Medical University of Gdansk, Gdansk, Poland

  • Brenda Schilke,

    Roles Investigation, Methodology, Writing – review & editing

    Affiliation Department of Biochemistry, University of Wisconsin-Madison, Madison, Wisconsin, United States of America

  • Rafal Dutkiewicz,

    Roles Data curation, Investigation, Methodology, Writing – review & editing

    Affiliation Intercollegiate Faculty of Biotechnology, University of Gdansk and Medical University of Gdansk, Gdansk, Poland

  • Marta A. Uzarska,

    Roles Investigation, Methodology, Writing – review & editing

    Affiliation Intercollegiate Faculty of Biotechnology, University of Gdansk and Medical University of Gdansk, Gdansk, Poland

  • Szymon J. Ciesielski,

    Roles Investigation, Methodology, Writing – review & editing

    Affiliation Department of Biochemistry, University of Wisconsin-Madison, Madison, Wisconsin, United States of America

  • Jacek Czub ,

    Roles Conceptualization, Data curation, Formal analysis, Funding acquisition, Project administration, Resources, Supervision, Validation, Visualization, Writing – original draft, Writing – review & editing

    jacek.czub@pg.edu.pl (JC); ecraig@wisc.edu (EAC); jaroslaw.marszalek@ug.edu.pl (JM)

    Affiliation Department of Physical Chemistry, Gdansk University of Technology, Gdansk, Poland

  • Elizabeth A. Craig ,

    Roles Conceptualization, Funding acquisition, Project administration, Resources, Supervision, Validation, Visualization, Writing – original draft, Writing – review & editing

    jacek.czub@pg.edu.pl (JC); ecraig@wisc.edu (EAC); jaroslaw.marszalek@ug.edu.pl (JM)

    Affiliation Department of Biochemistry, University of Wisconsin-Madison, Madison, Wisconsin, United States of America

  • Jaroslaw Marszalek

    Roles Conceptualization, Funding acquisition, Project administration, Resources, Supervision, Validation, Visualization, Writing – original draft, Writing – review & editing

    jacek.czub@pg.edu.pl (JC); ecraig@wisc.edu (EAC); jaroslaw.marszalek@ug.edu.pl (JM)

    Affiliations Intercollegiate Faculty of Biotechnology, University of Gdansk and Medical University of Gdansk, Gdansk, Poland, Department of Biochemistry, University of Wisconsin-Madison, Madison, Wisconsin, United States of America

Abstract

J-domain proteins (JDPs), obligatory Hsp70 cochaperones, play critical roles in protein homeostasis. They promote key allosteric transitions that stabilize Hsp70 interaction with substrate polypeptides upon hydrolysis of its bound ATP. Although a recent crystal structure revealed the physical mode of interaction between a J-domain and an Hsp70, the structural and dynamic consequences of J-domain action once bound and how Hsp70s discriminate among its multiple JDP partners remain enigmatic. We combined free energy simulations, biochemical assays and evolutionary analyses to address these issues. Our results indicate that the invariant aspartate of the J-domain perturbs a conserved intramolecular Hsp70 network of contacts that crosses domains. This perturbation leads to destabilization of the domain-domain interface—thereby promoting the allosteric transition that triggers ATP hydrolysis. While this mechanistic step is driven by conserved residues, evolutionarily variable residues are key to initial JDP/Hsp70 recognition—via electrostatic interactions between oppositely charged surfaces. We speculate that these variable residues allow an Hsp70 to discriminate amongst JDP partners, as many of them have coevolved. Together, our data points to a two-step mode of J-domain action, a recognition stage followed by a mechanistic stage.

Author summary

It is well appreciated that Hsp70-based systems are the most versatile among molecular chaperones—functioning in all cell types and in all subcellular compartments. Via cyclic binding to protein substrates, Hsp70s facilitate their folding, trafficking, degradation and ability to interact with other proteins. Hsp70 function, however, depends on transient interaction with J-domain protein cochaperones that not only deliver substrates, but also activate the structural changes needed for efficient Hsp70 binding to substrate. But how J-domain proteins mechanistically function to drive these changes and how an Hsp70 discriminates among multiple J-domain partners have remained challenging central questions. Here, by using a combination of computational, evolutionary and experimental approaches, we provide evidence for a two-step mechanism. The initial recognition step involves variable residues that allow fine tuning of both the specificity and strength of J-domain protein interaction with Hsp70. The second, that is the mechanistic step, involves conserved residues that act to disrupt a conserved network of intramolecular interactions within Hsp70, thus ensuring robust activation of the structural changes necessary for effective substrate binding. We suggest that our findings are likely applicable to most Hsp70 systems.

Introduction

By transiently binding many different polypeptide substrates, Hsp70 chaperones assist in diverse cellular processes, from de novo protein folding to protein trafficking to disassembly of protein complexes. J-domain protein (JDP) cochaperones are crucial players in the cycle of interaction of Hsp70 with substrate. They transiently interact with ATP-bound Hsp70 via their defining J-domain. This interaction, in coordination with substrate binding, facilitates stimulation of Hsp70's ATPase activity, which in turn drives the large scale conformational changes that stabilize substrate interaction [14].

The structure of both J-domains and Hsp70s are conserved. J-domains are composed of four α-helices. Helices II and III form a finger like structure. The helix II/III connecting loop includes a conserved histidine, proline, aspartate tripeptide (HPD). These invariant residues are critical for J-domain stimulation of Hsp70's ATPase activity. Hsp70s are more complex in structure, containing two large domains—a nucleotide binding domain (NBD) and a substrate binding domain (SBD)—connected by a “linker” segment. When ATP is bound, the domains are docked such that the substrate binding site in the β subdomain of the SBD (SBDβ) is easily accessible [5, 6]. Upon ATP hydrolysis, the domains disengage, and the substrate is “trapped” as the α subdomain of the SBD (SBDα) closes over the SBDβ substrate binding site [79].

The heart of the allosteric control mechanism behind these Hsp70 conformational transitions is an interaction network at the interdomain interface [1012]. A critical segment of this interface is formed by interactions at the base of the NBD with the linker and SBDβ (called NBD/SBDβ,linker throughout) [5, 6, 13]. The linker plays an important role, as its interaction with the NBD is key to stimulation of ATPase activity [1417]. A transient allosterically active intermediate poised for hydrolysis of ATP, often referred to as the “allosterically active state”, has been observed in the presence of a peptide substrate—NBD/SBD contacts are largely absent, but the linker remains bound to the NBD [10, 18].

It has been clear for years that J-domain interaction is essential for stimulation of Hsp70’s ATPase activity [19]. Particularly relevant to this report, early work on Escherichia coli JDP DnaJ and Hsp70 DnaK established an arginine at the Hsp70 NBD/SBDβ,linker interface, R167 of the NBD, as an important for J-domain function [14, 20, 21]. More recently a crystal structure of the J-domain of DnaJ in complex with ATP-bound DnaK was obtained [22]. This normally transient interaction was trapped by covalently linking the J-domain and Hsp70 in cis, and by crosslinking Hsp70’s NBD and SBD to stabilize intramolecular interactions between domains. Consistent with the results of a wealth of mutational and biochemical analyses, as well as earlier molecular dynamics simulations [23], the J-domain is bound at the NBD/SBDβ,linker interface. Residues of helices II and III and the connecting HPD loop form polar and hydrophobic interactions with the NBD and SBDβ.

Despite the accumulated information about J-domain/Hsp70 interactions, important questions remain. Although the end result of ATP hydrolysis is destabilization of Hsp70 domain-domain interactions, the initial structural consequences of J-domain binding are not known. Furthermore, there is no understanding of what determines the specificity of J-domain/Hsp70 interactions and how Hsp70s discriminate amongst their multiple JDP partners. To address these issues, we investigated two JDP/Hsp70 pairs: DnaJ/DnaK of E. coli and the mitochondrial Hsc20/Ssq1 of Saccharomyces cerevisiae [24, 25], which specializes in the biogenesis of proteins containing iron-sulfur clusters.

Results of our analyses point to two previously unappreciated facets of J-domain/Hsp70 interactions, leading us to propose a two-step model of J-domain action–recognition, followed by mechanistic action. Molecular simulations point to an active role of the invariant aspartate of the J-domain HPD in the mechanistic step–perturbing the intramolecular interaction network formed by the previously identified, key NBD arginine. This perturbation promotes allostery-related conformational changes of Hsp70. Evolutionary analyses indicate that J-domain/Hsp70 recognition, is dependent on variable, coevolving residues, even though the general mode of interaction of J-domains with Hsp70s is the same. This variability provides a plausible explanation for how Hsp70s are able to discriminate amongst their multiple JDP partners, thus enabling evolution of complex JDP/Hsp70 interaction networks.

Results

Structural model of the Hsc20-Ssq1 complex

To obtain a structural model of the JDP-Hsp70 complex formed by Hsc20 and Ssq1, we combined protein-protein docking with all-atom molecular dynamics (MD) simulations (see Methods and S1 Fig for details). This approach revealed a dominant bound state that accounted for 56% of the bound population (Fig 1A, S2 Fig). No other bound state accounted for more than 6% of the ensemble; these minor states rapidly interconverted and had a clear tendency to converge to the most populated bound state (S2 Fig). The binding mode in this dominant state closely resembles the recently published X-ray structure of DnaJ’s J-domain in complex with DnaK (DnaJJD-DnaK) [22]. Helices II and III are oriented very similarly in the two complexes—with a backbone root mean square deviation (RMSD) of 4.17 Å (Fig 1B). More specifically, the J-domains bind at the NBD/SBDβ,linker interface such that helix II predominantly interacts with the β-sheet region of NBD subdomain IIa and helix III with SBDβ (Fig 1). The HPD residues are positioned nearly identically in the two complexes, towards R167/R207 of NBD subdomain Ia in DnaK/Ssq1 (referred to as RNBD throughout), the arginine that previous analyses of DnaK established as important for allosteric transitions [14, 20, 21]. The predicted model of the Hsc20-Ssq1 complex was additionally validated using biochemical experiments, which are described later in the text.

thumbnail
Fig 1. Binding mode of J-domain-Hsp70 complexes.

(A) (top left) The dominant state of the Hsc20-Ssq1 complex obtained by protein-protein docking followed by MD refinement. Structure of S. cerevisiae Hsc20 (PDB ID 3UO3) was used along with a homology model of ATP-bound Ssq1, based on the X-ray structure of ATP-bound DnaK (PDB ID 4B9Q). The Ssq1 model had the expected architecture of an ATP-bound Hsp70: Nucleotide binding domain (NBD, grey) interacting with α and β subdomains of the substrate binding domain (SBD, brown) and the linker (yellow), in a groove between NBD subdomains Ia (light grey) and IIa (dark grey). To model the Hsc20-Ssq1 complex, conformational fluctuations of Hsc20 and ATP-bound Ssq1 were first characterized using MD simulations and clustering (S1 Fig). The obtained distinct conformers were then combinatorially docked. Based on the energetic and geometric criteria described in Methods, a representative set of 33 of the resulting complexes was selected for further refinement with MD simulations, and then clustered according to the positioning of the J-domain with respect to Ssq1 (S2 Fig). The model shown represents the dominant state with a population of 56%: J-domain of Hsc20 (cyan) interacts with Ssq1 at the NBD/SBDβ, linker interface. HPD sequence (red) of the J-domain is positioned near the conserved R207 (blue) of the NBD. (top middle) View of the Hsc20-Ssq1 complex to illustrate the orientation of the J-domain relative to the NBD subdomains (IA, IB, IIA, IIB s) which are marked, only the J-domain of Hsc20 is shown for clarity. (top right) Close-up showing J-domain interaction with an orientation similar to that at left. (bottom left) X-ray structure of the J-domain of DnaJ (DnaJJD) in complex with DnaK (DnaJJD-DnaK) (PDB ID 5NRO). (B) Structural alignment of the predicted Hsc20-Ssq1 complex and DnaJJD-DnaK crystal structure (the backbone RMSD of helices II and III of the J-domains after aligning NBDs is 4.17 Å). The J-domain, NBD, linker and SBD of aligned proteins are represented in cyan/magenta, grey/pink, yellow/green and light/dark orange for Hsc20-Ssq1 and DnaJJD-DnaK, respectively.

https://doi.org/10.1371/journal.pcbi.1007913.g001

DHPD interferes in Hsp70 interdomain interactions by perturbing critical RNBD contacts

To address the question of the mechanism of J-domain action, we employed MD simulations that, unlike most other methods, allow studying of transient residue-residue interactions [26]. We first analyzed Hsp70 alone, starting with the Ssq1 structural model (see Methods for details)—focusing on intramolecular interactions involving the RNBD because of its demonstrated importance in NBD/SBD communication [14, 20, 21] and its proximity to the HPD in the J-domain-Hsp70 complexes. Our initial analysis of MD trajectories revealed RNBD to be capable of forming two ion pairs, one with D517 of the SBDβ (referred to as DSBD throughout) and one with D429 of the interdomain linker (referred to as DLk throughout). To analyze these two contacts of RNBD more closely, we constructed the free energy map governing their formation (Fig 2A), by applying metadynamics, an enhanced sampling method suitable for simulating complex conformational changes [27, 28]. As relevant coordinates, we used the distances between the centers of mass of the guanidine group of RNBD and carboxylic groups of its interaction partners DSBD and DLk. Using the distance threshold for an ion pair of 0.5 nm, we found the probability of the two contacts occurring at the same time to be 82% (Fig 2A, S1 Table). In the DnaK crystal structure [5] the homologous DSBD- RNBD pair is present, but the homologous DLk-RNBD pair is not. To address this apparent discrepancy, we performed metadynamics simulations of the DnaK crystal structure. Analysis of the free energy map revealed that DnaK RNBD interacts simultaneously with DSBD and DLk 83% of the time. Taken together, our metadynamics simulations indicate that the RNBD double interaction can occur frequently in both Ssq1 and DnaK, despite the fact that it was not captured in the static DnaK crystal structure.

thumbnail
Fig 2. DHPD of J-domain displaces interdomain contacts formed by RNBD with the DLk and DSBD.

(a) Free energy landscape governing the formation of Ssq1/DnaK RNBD interdomain contacts with DLk and DSBD (top) Schematic depicting of the three Hsp70 residues of the RNBD interaction network; RNBD of subdomain Ia (blue dot) can interact with DLK of the linker (yellow dot) and DSBD of the SBDβ (brown dot). The four possible interaction states of RNBD are illustrated at right with the color of the rectangles serving as a key for the free energy landscapes at bottom. NBD subdomains Ia (light grey) and IIa (dark grey), SBD subdomains α and β (brown), linker (yellow). RNBD = R207/R167 in Ssq1/DnaK; DLk = D429/D393 in Ssq1/DnaK; DSBD = D517/D481 in Ssq1/DnaK. (bottom) Free energy landscapes governing the formation of these interaction states for Ssq1 and DnaK alone (left panels) and for Hsc20-Ssq1 and DnaJJD-DnaK complexes (right panels) with distances between RNBD—DLk and between RNBD- DSBD used as relevant coordinates. Integrated populations of these interaction states are shown by pie charts (see S1 Table for data with standard errors). (b) Free energy profile of the transition of RNBD from intramolecular interactions with DLk and DSBD residues to intermolecular contact with DHPD of the J-domain HPD loop when kept near NBD subdomain Ia. Distance difference Δd = d1—d2 is used as a relevant coordinate. d1 and d2 are the distances between RNBD and DLk/DSBD of Hsp70 (d1) or DHPD of the J-domain (d2). Conformations in which RNBD interacts with DHPD at the expense of interaction with DLk and/or DSBD (representative structures on the right, corresponding to the free energy minimum at 0.4–0.5 nm) are energetically more favorable than conformations with RNBD interacting simultaneously with DLk and DSBD (representative structures on the left).

https://doi.org/10.1371/journal.pcbi.1007913.g002

To assess the effect of J-domain binding to Hsp70 on this conserved arginine-centered interaction network, we carried out the same set of calculations for Ssq1 and DnaK in complex with Hsc20 and DnaJJD, respectively (Fig 2A, S1 Table). For both Hsc20-Ssq1 and DnaJJD-DnaK, the fraction of the population of Hsp70 molecules in which RNBD contacts both DSBD and DLK was reduced compared to that of Hsp70 alone: Ssq1, from 82 to 62%; DnaK, from 83 to 51%. We hypothesized that DHPD competes with DSBD and DLk for interactions with RNBD and therefore is responsible for the observed partial loss of its intramolecular contacts. We also used metadynamics trajectories to calculate all RNBD contacts both in the presence and in the absence of the J-domain (S3 Fig); beside the aspartates discussed above, only interaction with the Ala residues next to DSBD, was reduced in the presence of the J-domain. To determine the likelihood of such an exchange of interaction partners, we performed metadynamics simulations to compute the free energy profiles for transition of RNBD from intramolecular interaction with DSBD and DLk to interaction with DHPD (Fig 2B). In these simulations, the HPD loop was kept near NBD subdomain Ia in the orientation predisposing it for the interaction with RNBD (see Fig 1A). The well-defined free energy minima in the resulting profiles indicate that RNBD shows preference (1.5 and 4 kcal/mol for Ssq1 and DnaK, respectively) for interaction with DHPD, rather than simultaneous interaction with the aspartates on the SBD and the linker, DSBD and DLk. The observed differences in the depth of the free energy minima between Ssq1 and DnaK could be attributed to differences in the amino acid sequences or to the fact that Hsc20 was modeled as a whole instead of only its J-domain. Finally, an example of spontaneous transition of RNBD to the intermolecular contact with DHPD, captured in our unbiased equilibrium MD simulation of the Hsc20-Ssq1 complex, is shown in S1 Movie. Taken together, our data suggest that DHPD perturbs the RNBD -centered interaction network at the NBD/SBDβ,linker interface.

J-domain binding initiates SBDβ disengagement from the NBD

We next asked how perturbation of intramolecular interactions of RNBD with DSBD and DLk affects the NBD/SBDβ,linker interface. First, we considered the effect of substitution of RNBD, such that interactions with aspartates are prohibited. An RNBD->D DnaK substitution variant has been shown to be defective in interdomain communication [14]. Its intrinsic ATPase activity is somewhat higher than that of wild-type (WT) DnaK and not effectively stimulated in the simultaneous presence of its JDP partner DnaJ and protein substrate. We therefore asked whether substitution of RNBD->D in Ssq1 had a similar effect. We found the basal ATPase activity of Ssq1 RNBD->D to be 3-fold higher than that of WT Ssq1, and not efficiently stimulated in the presence of Hsc20 and protein substrate, reaching only 30% of the activity observed for WT Ssq1 (S4 Fig). Thus, substitution of the conserved RNBD has significant effects on interdomain communication in both DnaK and Ssq1.

Because of these effects on interdomain communication, we used MD simulations to probe how weakening of the cross-domain interactions of RNBD (i.e. with DSBD and DLk) affects the compactness of the NBD/SBDβ,linker interface. Umbrella sampling was used because it allows testing of conformations not accessible by conventional MD simulations. For these simulations a center of mass distance between SBDβ and its NBD interacting region (SBDβ-NBD distance) was used as a relevant coordinate. To accelerate the free energy convergence and to allow exhaustive sampling of interdomain contacts these simulations were performed using Ssq1 and DnaK with a truncated SBDα subdomains.

Free energy profiles were computed (S5 Fig), and then converted into probability distributions (Fig 3A). These distributions show that Ssq1 and DnaK both adopt conformations with longer SBDβ-NBD distances in the presence, than in the absence of a J-domain, indicative of a transition to a less compact interdomain interface. Furthermore, the probability distributions of the RNBD->D variants (Fig 3A) resemble those of Hsp70s in the presence of the J-domain in that that the mean distance values are very similar for both distributions, supporting our hypothesis that disruption of the arginine centered interaction network affects compactness of the interdomain interface.

thumbnail
Fig 3. J-domain binding displaces SBDβ from association with NBD.

(a) Structural representation of the center of mass distance, d, between SBDβ (orange) and its NBD interface (green). Probability distributions of the interdomain distance d for Ssq1 and DnaK; Ssq1 and DnaK alone (top panels, in grey), Hsc20-Ssq1 and DnaJJD-DnaK complexes (center panels, in cyan), RNBD substitution variants Ssq1(R207D) and DnaK(R167D) (bottom panels, in red). Dashed lines represent an average d value for each distribution. (b) The number of contacts between NBD and linker and between NBD and SBDβ, with increasing distance d; Ssq1 and DnaK (grey dots); Hsc20-Ssq1 and DnaJJD-DnaK (cyan dots); Ssq1(R207D) and DnaK(R167D) (red dots). Contacts were defined as the number of non-hydrogen atoms within a 0.5 nm cutoff distance across the binding interface.

https://doi.org/10.1371/journal.pcbi.1007913.g003

To assess what structural changes took place within the NBD/SBDβ, linker interface we used umbrella sampling trajectories to determine the number of interdomain contacts (Fig 3B). We analyzed interactions between the NBD and SBDβ and between the NBD and the linker separately, because of their opposing effects on ATPase activity (i.e. ATPase activity requires linker association and SBD dissociation). In all our analyses, for both Ssq1 and DnaK, the number of the NBD/linker contacts remains virtually the same regardless of the interdomain distance, while the number of NBD/SBDβ contacts decreases steadily with increasing distance (Fig 3B). Taken together, these data indicate that in the presence of either the J-domain or the RNBD->D substitution (i.e., when the interdomain distance increases), both DnaK and Ssq1 adopt conformations with markedly reduced number of NBD/SBDβ contacts. In particular, upon J-domain binding, the populations of DnaK and Ssq1 that retain at least 70% of interdomain contacts are reduced from 88% to 50% and from 91% to 45%, respectively. A similar reduction is observed for the RNBD->D substitution in the absence of J-domain, indicating that perturbation of the arginine-centered interaction network promotes partial disengagement of SBDβ. These results also indicate that J-domain binding perturbs the RNBD intramolecular interactions, which in turn leads to partial disengagement of SBDβ, while the interdomain linker remains tightly associated with the NBD. Such an arrangement favors ATP hydrolysis, which triggers further conformational changes [10, 11, 14, 15].

Conservation of residues involved in the mechanistic step of J-domain action

The RNBD-centered interaction network is the same in DnaK and Ssq1, in that RNBD interacts with homologous aspartates in the linker and SBDβ. To investigate the prevalence of these residues, we analyzed Hsp70 sequences from a large number of bacterial and eukaryotic species. RNBD and its interaction partner DLk is invariant in our data set (Table 1). The conservation of its SBDβ interaction partner DSBD is somewhat more complex. Aspartate is only present at this position in DnaKs of bacteria closely related to E. coli and in mitochondrial Hsp70s of fungi; asparagine is present at this position in most Hsp70s (Table 1, S6 Fig). That asparagine can also form a hydrogen bond with arginine [29], thereby participating in the interaction network, supports the idea that the arginine-centered network is conserved across Hsp70s. Thus, it is likely the mechanism of action of J-domains, that is perturbation of this network, is conserved as well.

thumbnail
Table 1. Evolutionary conservation of RNBD-centered interaction network.

https://doi.org/10.1371/journal.pcbi.1007913.t001

Importance of electrostatic interactions at the Hsc20/Ssq1 binding interface

Having found similarities of the RNBD network and the overall mode of interaction of DnaJJD-DnaK and Hsc20-Ssq1, we decided to look more closely at the J-domain/Hsp70 binding interfaces. DnaJJD-DnaK has been studied extensively [22]; both polar and hydrophobic contacts were identified based on the DnaJJD-DnaK crystal structure and the functional importance of key residues was verified experimentally. However, no information is available for Hsc20-Ssq1. We therefore expanded our computational analysis of the Hsc20-Ssq1 complex. We used MD simulation trajectories to search for all specific (polar and hydrophobic) residue-residue contacts, and computed the relative energetic contributions of individual residues to complex stability (S7 Fig, S8 Fig). The binding interface is formed by 17 residues of the J-domain contacting 18 residues of Ssq1 at a cutoff distance of 0.5 nm. The interface, which occupies 6.7 ± 0.2 nm2, is mostly hydrophilic, with scattered hydrophobic interactions contributing 12.5% of the surface (0.8 ± 0.1 nm2). The predominant contributors are eight charged residues (Fig 4A). Three positively charged residues of Hsc20 helix II (R37JD, K38JD, R41JD) that form a network of electrostatic interactions with four negatively charged residues of the Ssq1 NBD (D246, E248, D249, E253) and K70JD of helix III that promotes complex formation, though its binding partner(s) was not identified in our analyses.

thumbnail
Fig 4. Electrostatic interface of Hsc20 J-domain and Ssq1.

(a) Network of key polar interactions at the binding interface between the J-domain of Hsc20 (cyan) and Ssq1 subdomain IIa (grey), linker (yellow) and SBDβ (brown). The most stable ion pairs and hydrogen bonds identified in MD simulations are indicated with dashed lines. Relative contributions of individual residues to the Hsc20-Ssq1 complex stability are indicated by color scale; blue for J-domain residues and red for Ssq1 residues. Note: K70JD does not form obvious contacts but contributes strongly to the complex stability. (b) Site specific crosslinking of Hsc20 to Ssq1. Purified Hsc20 Q46Bpa (left) and M51Bpa (right) variants were incubated with purified Ssq1 wild-type (WT) or alanine substitution variants (D246A, E253A or E248A, D249A) in the presence of ATP. After UV irradiation (+), or as a control no irradiation (-), reaction mixtures were separated by SDS-PAGE. Migrations of size standard, in kDa, are indicated; Q:J* and Q:J** indicate positions of the Ssq1-Hsc20 crosslinks, which were identified using mass-spectroscopy (S9 Fig); we do not understand the nature of the more slowly migrating crosslink band (Q:J**) present in the M51Bpa variant sample. Quantification of three independent Bpa crosslinking experiments; statistically significant differences (t-student test, p<0.001) are indicated (***). (c) Stimulation of Ssq1 ATPase activity by Hsc20 was measured in the presence of 0.5 μM of Ssq1, 500 μM of PVK-peptide substrate, 0.5 μM of nucleotide exchange factor Mge1 and indicated concentrations of Hsc20. (top) Hsc20 WT and alanine substitution variants. (bottom) Ssq1 WT and alanine substitution variants. Error bars indicate standard deviation for three independent measurements.

https://doi.org/10.1371/journal.pcbi.1007913.g004

We carried out two types of biochemical experiments to verify the functional importance of key residues that based on our computational analysis form the Hsc20/Ssq1 binding interface. First, we asked if alanine substitution of these Ssq1 residues disturbs the Hsc20/Ssq1 interaction. To enable detection of the inherently transient J-domain/Hsp70 interaction, we used Hsc20 having a site-specific p-benzoyl-l-phenylalanine (Bpa) crosslinker [30] replacing two residues, one upstream (Q46) and one downstream (M51) of the HPD (Fig 4B, S9 Fig). The efficiency of Hsc20 crosslinking to the Ssq1 variants was strongly reduced compared to that of WT Ssq1. As a control for crosslinking specificity we replaced a residue outside the binding interface with Bpa (S34Bpa); no crosslink band(s) was detected (S10 Fig). Second, we measured Hsc20’s ability to stimulate the ATPase activity of Ssq1. We tested four Hsc20 and four Ssq1 variants, each having a different single alanine substitution. All variants had reduced activity (Fig 4C). We conclude that the key interfacial residues of Hsc20 and Ssq1 identified by the MD simulations are indeed functionally important for the Hsc20/Ssq1 interaction.

Hsc20 binding to Ssq1 is driven by long-range electrostatic interactions

The functional importance of charged residues in Hsc20/Ssq1 interaction led us to consider whether the initial encounter of the two proteins is electrostatically driven—an idea furthered by the presence of a patch of uniformly positive electrostatic potential around helix II of the J-domain and a complementary patch of negative potential around subdomain IIa of the Ssq1 NBD, the site of J-domain interaction (S11 Fig). MD simulations in which Hsc20 and Ssq1 were initially separated with their charged patches facing each other, revealed that Hsc20 spontaneously interacted with Ssq1 within tens of ns and remained stably bound for the rest of the stimulation (S2 Movie). Notably, the spontaneously formed complexes are very similar to the dominant bound mode identified in our systematic search (S11 Fig).

To examine the driving forces behind spontaneous Hsc20/Ssq1 interaction, we computed the underlying free energy profile using umbrella sampling. The resulting profile (S11 Fig, shows a pronounced bound-state energy minimum (at distances < 2 nm), indicative of site-specific binding. The non-zero slope of the WT free energy profile up to ~3.5 nm demonstrates that Hsc20 is attracted to Ssq1 by long-range electrostatic driving forces. Consistent with the above analysis of the binding interface and biochemical data (Fig 4), the free energy minimum, and hence binding affinity, almost completely disappears in the case of the Hsc20 variant having helix II alanine substitutions R37A,R41A (S11 Fig).

Residues forming the JDP/Hsp70 interface are variable and coevolving

The two divergent J-domain/Hsp70 partnerships analyzed here show remarkable structural similarity in their overall mode of interaction—J-domain helix II with the β-sheet of NBD subdomain IIa and J-domain helix III with SBDβ (Fig 1). However, when we aligned amino acid sequences of these regions it became apparent that the residues forming the two J-domain/Hsp70 binding interfaces are not evolutionary conserved (Fig 5A). This is particularly evident for the J-domains—only 12 of the 17 positions forming the Hsc20/Ssq1 interface are shared with the DnaJJD/DnaK interface, and of these 12, only 3 are occupied by the same amino acids in the two interfaces. A broader comparison of Hsc20, Ssq1, DnaJ, and DnaK with their orthologs from a wide range of species belonging to fungi, animals and plants (S12 Fig, S13 Fig) revealed that the amino acids occupying the J-domain/Hsp70 interfacial positions are indeed variable. The J-domain side of these interfaces is particularly so. Amino acid sequence variability is evident for J-domains of Hsc20 and DnaJ among relatively closely related species, Ascomycota and Bacteria, respectively, suggesting that even closely related JDP/Hsp70 pairs engage different sets of residues to form the binding interface (Fig 5B). One trend is evident however. Despite the overall sequence variability, the J-domain side of the interface is enriched in positively charged amino acids (Fig 5C)—consistent with our results that complementary electrostatic interactions are crucial for stability of the Hsc20-Ssq1 complex. The Hsp70 side of the binding interface is less variable than the J-domain side; 8 out of the 16 positions shared between the Ssq1/DnaK interfaces are occupied by identical amino acids.

thumbnail
Fig 5. Sequence variability of J-domain/Hsp70 binding interface.

(a) Pairwise sequence alignments of the J-domain/Hsp70 interfacial regions of the J-domains Hsc20 and DnaJ’s (top) and the Hsp70s Ssq1 and DnaK (bottom). Residues of each protein that are in contact across the Hsc20/Ssq1 and DnaJJD/DnaK binding interfaces are in green; HPD motifs in purple. (b) Sequence variability of positions that form the J-domain/Hsp70 binding interfaces (those positions shown in green in panel a) for orthologs of Hsc20/Ssq1 from 67 Ascomycota species (top) and DnaJ/DnaK from 117 bacteria species (bottom). Coevolving positions with p>0.01 statistical support (S2 Table, S3 Table, S16 Fig and S17 Fig) are connected by dashed lines. Sequence logos represent the amino acid frequency of each interfacial position in the orthologs examined, with positively charged (blue), negatively charged (red) and uncharged (grey) residues. Position numbering is that for S. cerevisiae (top) and E. coli (bottom) sequences. (c) Fraction of interfacial Hsc20/Ssq1 and DnaJJD/DnaK residues that are positively and negatively charged. Whiskers indicate standard deviation.

https://doi.org/10.1371/journal.pcbi.1007913.g005

The data discussed in the previous paragraph presents an apparent conundrum–sequence variability of J-domain/Hsp70 interfacial residues, yet maintenance of the functional interaction of J-domains with Hsp70 partners. Coevolution of residues across the binding interface, such that amino acid substitutions on one side are compensated by complementing substitutions on the other [31], could be an explanation.

Therefore, we tested whether J-domain/Hsp70 interfacial residues are coevolving. We analyzed the Ascomycota and bacterial data sets discussed above (Fig 5B), applying the probabilistic Coev model [32] of sequence evolution that allowed us to discriminate between coevolving and independently evolving positions, while incorporating the underlying phylogenetic relationships among analyzed sequences (S14 Fig, S15 Fig). We found that some of the variable positions are indeed coevolving—16 for the Hsc20/Hsp70 interface; 20 for DnaJ/DnaK (Fig 5B). As a control, we tested the same datasets using a model of 'independent' sequence evolution. The Coev model fits our data significantly better than the independent model (S2 Table, S3 Table, S16 Fig and S17 Fig), supporting our hypothesis that coevolution of residues across the binding interfaces could explain how J-domain/Hsp70 functional interactions are maintained despite sequence diversity.

Discussion

Based on the results presented we propose a mechanism of J-domain function that reconciles two apparently contradictory features: an invariant dependence on the HPD signature motif, yet an ability to recognize a specific Hsp70 partner. Our data supports a model in which the DHPD of the J-domain plays an active role in perturbing the intramolecular interactions of the Hsp70 RNBD with conserved residues of the linker and SBDβ, destabilizing NBD/SBD interactions and thus promoting allosteric transition leading to ATP hydrolysis. On the other hand, evolutionarily variable J-domain residues are responsible for the specificity of exclusive JDP/Hsp70 interactions, as well as fine-tuning prioritization of JDP/Hsp70 interactions in the cases where an Hsp70 has multiple JDP partners, each of them guiding Hsp70 to different cellular functions.

Mechanistic step

Based on detection of an energetically favorable interaction between invariant DHPD and invariant RNBD, we propose that the DHPD—RNBD interaction displaces those of RNBD with residues of the interdomain linker (DLk) and SBDβ subdomain (DSBD) of Hsp70. This displacement promotes destabilization of the interdomain interface, such that the SBDβ partially disengages from the NBD. These results are consistent with experimental data showing less protection against hydrogen exchange of this region when residues in the NBD/SBDβ interface are mutated [33] It is important to note that the interdomain linker remains stably bound to the NBD upon J-domain induced disengagement of the SBDβ. This structural arrangement has been termed the "allosterically active state" [10] because linker engagement with, but SBDβ disengagement from, the NBD is necessary for ATPase activation that triggers further Hsp70 conformational changes required for stabilization of substrate binding [11, 1417, 34, 35]. The data presented here also allows a more unified picture of how Hsp70 reaches the allosterically active state to effectively trap substrate. It was previously shown that the allosterically active state is populated in the presence of ATP and substrate [10, 11, 36]. Substrate binding induces an allosteric signal, which is propagated through a cascade of interactions from the SBDβ substrate binding site to the NBD/SBDβ,linker interface. Our data indicates that the allosterically active state it is also populated in the presence of ATP and J-domain. In this case, the signal is propagated by perturbation of the RNBD centered interaction network. These perturbations could also modify the interaction mode of ATP bound to Hsp70, thus affecting its rate of hydrolysis—a hypothesis that requires further studies to verify. Overall, we postulate that the substrate and J-domain dependent allosteric pathways converge and that this convergence allows synchronization—only upon simultaneous interaction of substrate and J-domain is ATPase activity effectively stimulated, driving transition to the undocked state and substrate capture.

This coordinated mechanism also impacts the larger picture of Hsp70 allosteric dynamics. All Hsp70s studied share major characteristic transient conformational shifts–for example, undergoing repeated transitions involving undocking from the NBD of both SBDβ and SBDα. However, recent single molecule experiments indicate that the frequency of such transitions differs among members of the Hsp70 family [79, 18, 3739]. Interestingly, eukaryotic Hsp70s such as BiP of the endoplasmic reticulum [8, 38] or Hsp70 and Hsc70 [39] of the cytosol, display conformational dynamics markedly different from DnaK and mitochondrial Hsp70. They have DSBD replaced by asparagine, pointing to the importance of the RNBD centered interactions for the allostery of Hsp70s.

Though no model for the mechanistic role of the DHPD- RNBD interaction per se has been previously put forth, the functional connection between these residues is well supported by genetic and biochemical data. Alteration of either RNBD or residues interacting with it (DSBD and DLk) [10, 11, 14, 17, 36, 40] disrupt allosteric communication between domains, resulting in a somewhat increased basal activity and ineffective stimulation by its JDP partner, DnaJ in the case of DnaK [14] and Hsc20 in the case of Ssc1 (S4 Fig). Substitution of RNBD(R167H) of DnaK suppresses nonfunctional DHPD(D35N) variant of DnaJ [20]. However, we do not mean to say that this role of the DHPD is the only mechanistic effect of J-domain binding to Hsp70. For example, the invariant histidine and/or proline of the HPD may play specific roles, and, as previously suggested [41], conformational strain induced in the J-domain upon binding to Hsp70 might be important as well. However, there is no question that RNBD is key and could be considered as an allosteric perturbation site for the J-domain action.

Recognition step

At first glance, the sequence diversity of the interface of J-domain-Hsp70 complexes, but maintenance of overall positioning, and thus, function, seems surprising. It raises the question of how this diversity evolved without loss of function. Our analyses suggest that functional J-domain/Hsp70 interactions have been maintained over time via coevolution of residues that form the binding interface. The nature of the Hsc20/Ssq1 interface is instructive in this regard. Each key residue on one side of the large interface contacts more than one residue on the opposite side. Such a multifaceted interaction network is highly flexible. Over time residues can be added to or deleted from the interface without ever losing the ability to maintain functional interaction between the coevolving partners, enabling modulation of both the specificity and the strength of the interaction. Coevolution could also explain how a JDP might have been able to switch Hsp70 partnership during evolution, such as occurred in yeast mitochondria—switching of Hsc20 from interacting with a general multifunctional Hsp70 (Ssc1) to a specialized Hsp70 (Ssq1) [42]. A JDP that interacts with a given Hsp70 partner could gain the ability to interact with a “new” Hsp70 (i.e. product of gene duplication), upon emergence of mutation(s) that facilitate a new interaction, without affecting the initial one, at least during the period of transition [43].

What effect might this diversity have had on the functionality of Hsp70 chaperone systems more generally? We suggest that it has allowed evolution of complex Hsp70/JDP networks, both in situations in which a single type of Hsp70 interacts with many different JDPs and in which more than one Hsp70 coexist in a cellular compartment. In the first situation, replacement could well modulate the strength of a JDP/Hsp70 interaction, explaining why some JDPs have "easier access" to their partner Hsp70 than others. Here, electrostatic interactions could play a particularly important role. A JDP with a highly charged J-domain binding face would be able to recognize its Hsp70 partner at a longer distance, giving it priority over others in forming a productive interaction. In the second situation, the diversity serves to "insulate" Hsp70 networks, each Hsp70 interacting effectively only with its own subset of JDPs. For example, it could explain how yeast cytosolic Hsp70s Ssa and Ssb, which share common ancestry, evolved into independent networks, each with its own JDP partners [44, 45].

In the two examples studied here, the predominant J-domain interactions are via helix II. However, their binding interfaces differ—the specialized Hsc20/Ssq1 interface involving mostly electrostatic interactions and the DnaJJD/DnaK interface being a mosaic of electrostatic and hydrophobic interactions [22, 41, 46]. However, across JDP/Hsp70s partnerships the interfaces may be even more varied. Interestingly, in two cases, Helix III of the J-domain was found to be very important for interaction with its Hsp70 partner—auxilin [35, 47] and polyomavirus T Antigen [48]. As only a few J-domain/Hsp70 interfaces have been studied, the breath of sequence variability, that supports proper positioning of the HPD sequence, the key to functional partnership, may be larger still.

Concluding remarks

Overall our data suggests a two-step process of J-domain function: initial binding that forms similar poses via highly variable residue-residue contacts, followed by a mechanistic step in which the invariant HPD perturbs the Hsp70 NBD-SBDβ,linker interface to foster a key allosteric transition. This separation could have future importance beyond a basic understanding of molecular chaperone function. It has long been appreciated that Hsp70 systems impinge on many biological functions and thus, not surprisingly, connected to many pathologic conditions raging from cancer, metabolic diseases to aging [49, 50]. As it has become appreciated that many Hsp70s partner with multiple JDPs [1, 45], the idea of targeting a JDP, rather than Hsp70 itself, has emerged [51]. Our findings showing a clear separation between recognition and mechanistic function provide a possible avenue to approach such intervention–modulation of the recognition step without affecting the mechanistic step, thus, manipulating the specificity or strength of particular J-domain/Hsp70 interactions.

Methods

Molecular Dynamics (MD) Simulations

All simulations were performed using Gromacs 5 [52] with the Plumed 2.1 plugin [53]. If not stated otherwise, the CHARMM36 force field [54] was used for proteins, ions and Mg-ATP, and the TIP3P model was used for water. In each of the simulation boxes, the numbers of Na+ and Cl- ions were adjusted to 0.15 M. Temperature was kept at 310 K with the v-rescale algorithm [55] using a coupling constant of 0.1 ps. Pressure was kept at 1 bar using the Parrinello-Rahman algorithm [56] with a coupling time of 5 ps. Periodic boundary conditions were applied and the Particle Mesh Ewald summation [57] was used to calculate long-range electrostatic interactions with a cut-off radius of 1 nm and a Fourier grid spacing of 0.12 nm. Van der Waals interactions were calculated with Lennard-Jones potential with a cut-off radius of 1 nm. All bonds involving hydrogen were constrained using the LINCS algorithm. Leap-frog Verlet algorithm was used to integrate equations of motion with a time step of 2 fs.

Ssq1 homology model

Ssq1 model was built with I-TASSER [58] using the ATP-bound structure of DnaK (PDB ID 4B9Q) as a template with C-score 0.34. The Mg-ATP complex and structural water molecules were placed in the conserved ATP binding pocket of Ssq1 according to the DnaK structure (PDB ID 4B9Q). The obtained model of Ssq1 was refined by 3-μs-long MD simulation with position restraints applied to Mg-ATP and structural water during the first 1 μs to ensure a proper relaxation of the ATP-binding pocket.

Model of Hsc20-Ssq1 complex

The structures of Hsc20 (PDB ID 3UO3) and Ssq1 (homology model) used for molecular docking were first subjected to unbiased MD simulations to sample global protein structural dynamics. The Hsc20 simulation system was composed of a single protein in a 9.2 nm × 9.2 nm × 9.2 nm dodecahedron box filled with ~17000 water molecules, while the Ssq1 simulation system was composed of a single ATP-bound protein solvated with ~48,500 water molecules in a 13 nm × 13 nm × 13 nm dodecahedron box. Both systems were simulated using CHARMM36 [54] for 10 and 2 μs, respectively, and independently, using AMBER99SB-ILDN [59] for 5 and 2 μs, respectively, to test for force-field dependence. The obtained MD trajectories were separately clustered with g_cluster using a 0.35 nm RMSD cut off for Hsc20 and a 0.7 nm RMSD cut off for Ssq1 (to account for large conformational fluctuations of the sub-domain SBDα). Centroids of each of the 12 and 7 clusters identified for Hsc20 and Ssq1, respectively, were treated as representative structures of the proteins and paired together in 84 independent docking runs using ClusPro [60]. As a result, 8,973 models were obtained. For each of the 84 docking runs, we selected the four best scoring models in each of the following criteria: a) binding energy estimation, b) minimal distance between the J-domain helix II and the NBD, c) minimal distance between the J-domain helix III and the NBD, d) minimal distance between the J-domain and the NBD. Redundancy among 336 selected models was removed with g_cluster, which resulted in 33 distinct models of the Hsc20-Ssq1 complex. Next, each of these models was placed in a 16 nm × 16 nm × 16 nm dodecahedron box and solvated with ~85000 water molecules. The resulting 33 systems were energy-minimized and the protein side chains were relaxed with 50 ns-long MD simulations, during which the protein backbone atoms were restrained. Finally, the systems were simulated by unbiased MD for a minimum of 500 ns (total time for 33 simulated systems was 21 μs). The resulting trajectories were clustered based on heavy atoms of helices II and III of the Hsc20 J-domain after superimposing the Hsc20-Ssq1 complex structures using all Cα atoms of the conformationally rigid NBD domain. This was done using g_cluster with a single linkage method and a RMSD cutoff of 0.18 nm to obtain a fine characterization of the bound-state ensemble.

MD simulations of DnaK and the DnaJJD-DnaK structures

The systems were composed of either DnaK (PDB ID 4B9Q) or DnaK in complex with the J-domain of DnaJ (PDB ID 5NRO) placed in 13 nm × 13 nm × 13 nm dodecahedron box and solvated with approx. 52,000 water molecules. Crystal water molecules coordinated to Mg2+ ion in Mg-ATP were kept. Both systems were simulated for 1 μs and the final structures from these simulations were used as initial configurations for the subsequent free energy computations.

Energy landscape of intramolecular contacts with the conserved RNBD upon J-domain binding to Hsp70

The free energy landscapes describing the formation of intramolecular ion pairs by the critical RNBD (R207 in Ssq1 and R167 in DnaK) were computed with well-tempered multiple walkers 2D-metadynamics [61]. In order to examine the effect of the J-domain binding on the RNBD interactions, the simulation systems contained either Hsc20-Ssq1 or DnaJJD-DnaK complexes or individual proteins, Ssq1 or DnaK, respectively, placed in a 11.8 nm × 11.8 nm × 11.8 nm dodecahedron box filled with ~34,000 water molecules. To make the system significantly smaller and thus achieve better convergence of the free energy maps, the lid sub-domain, SBDα, was truncated by removing residues 570–657 of Ssq1 and 535–602 of DnaK. As two coordinates for 2D-metadynamics, we used the center of mass (COM) distances between the guanidine group of the RNBD residue and the carboxylic groups of the partner aspartate residues at the SBDβ (D517 in Ssq1 or D481 in DnaK) and at the interdomain linker (D429 in Ssq1 or D481 in DnaK). Biasing potential was added every 10 ps with the bias factor of 15, a height of 0.1 kJ and a width of 0.04 nm. The sampled range of both coordinates was restricted to a 0.3–1.9 nm interval by applying one-sided harmonic potentials with a force constant of 3000 kJ nm-2. For both systems, 10 independent 1-μs-long walker-simulations were carried out. The final free energy maps were averaged over 15 individual maps computed in 25-ns intervals from the last 350 ns of the simulations to average out free energy fluctuations [62].

The free energy profiles for the transition of the conserved arginine residue from intramolecular ion pairs with SBDβ (D517 for Ssq1, D481 for DnaK) and the linker (D429 for Ssq1, D481 for DnaK) to an intermolecular contact with the aspartate in the HPD (D50 for Hsc20, D35 for DnaJ) were computed using well-tempered multiple walkers metadynamics. The sampled coordinate was defined as the difference between the two distances: i) the COM distance between the arginine (RNBD) guanidine group and the carboxylic groups of the aspartate residues in the linker (DLk) and in the SBDβ (DSBD) and ii) the COM distance between the RNBD guanidine group and the carboxylic group of the DHPD (see Fig 2). To speed up convergence, both of these component distance coordinates were restricted to up to 1.5 nm by applying one-sided harmonic potentials with a force constant of 3000 kJ nm-2. The bias potential was incremented every 10 ps by adding Gaussian-shaped potentials with height and width of 0.06 kJ and of 0.04 nm, respectively, and a bias factor of 8. For both systems, 10 walkers of 400 ns each were simulated at the same time. The final free energy profiles were averaged over 8 individual profiles computed in 20-ns intervals from the last 140 ns to limit free energy fluctuations.

Conformational transition of Hsp70 upon J-domain binding

The free energy profiles describing the compactness of the NBD, linker and SBDβ interface in Ssq1 and DnaK in the presence and absence of J-domain were computed with replica-exchange umbrella sampling (REUS). The simulations of the protein complexes (Hsc20-Ssq1 and DnaJJD-DnaK) and free proteins (Ssq1 and DnaK) were performed in a 11.8 nm × 11.8 nm × 11.8 nm dodecahedron box filled with ~36,000 water molecules. The reaction coordinate was defined as the COM distance between the Cα atoms of SBDβ (residues 434–541 in Ssq1 and 395–505 in DnaK) and the Cα atoms of the NBD region at the interface with SBDβ (residues 112–143 and 186–262 in Ssq1, residues 75–104 and 146–226 in DnaK). To accelerate the free energy convergence, the lid sub-domain (SBDα) was truncated (residues 542–657 in Ssq1 and 506–602 in DnaK). For all systems, the initial configurations for REUS simulations were alternately selected from two independent 200 ns steered-MD simulations. In these simulations, the distance between SBDβ and NBD was gradually increased using a moving one-sided harmonic potential with a force constant of 2500 kJ nm-2, applied to either the COM distance or the minimal distance between the Cα atoms of SBDβ and the NBD interfacial residues. The reaction coordinate in the 2.5–2.9 nm range was divided into 5 equally spaced windows in which the system was restrained using a harmonic potential with a spring constant of 2500 kJ nm-2. The same approach was used to study the effect of the RNBD->D substitutions (Ssq1(R207D) and DnaK(R167D)) on the relative arrangement of NBD and SBDβ. To perform these simulations, the RNBD residue was replaced by aspartate in all initial structures. For each of the REUS windows at least 800 ns long simulation were performed. Exchanges between neighboring umbrella sampling windows were attempted every 2 ps. The first 100 ns of the trajectory was discarded and free energy profiles were computed using WHAM [63] with the Monte Carlo bootstrap method to estimate uncertainties of the free energy. The number of contacts between the domains (NBD/SBDβ or NBD/linker) is defined as the number of non-hydrogen atoms within a cutoff distance of 0.5 nm. To compute the number of contacts as a function of the COM distances between the domains, REUS trajectories were reweighted using a factor of exp((Ui(r)−Fi)/kT), where Ui(r) denotes the value of the biasing potential for a given frame and Fi denotes the free energy added to the i-th US window by the applied bias.

Hsc20-Ssq1 interaction analysis

All inter-protein interactions were analyzed using a 10.5 μs MD trajectory of the Hsc20-Ssq1 complex representing a dominant bound state. (i) Ion pairs between Hsc20 and Ssq1 were identified using a cut-off of 0.5 nm for the distance between the centers of mass of the interacting charged groups [64]. (ii) Hydrogen bonds were identified using a standard geometric criterion: a 0.35 nm cut-off for donor-acceptor distance and a 30° cut-off for hydrogen-donor-acceptor angle. (iii) The hydrophobic contact surface area was estimated using the g_sas tool from the GROMACS package [52]. The hydrophobic surface buried in the Hsc20/Ssq1 interface was computed as the difference between the surface of hydrophobic residues exposed to water in the dissociated and bound state. (iv) Relative contributions to the Hsc20-Ssq1 binding free energy of individual residues were calculated with g_mmpbsa [65] using a single trajectory approach, non-linear Poisson-Boltzmann solver and the CHARMM charges and atomic radii.

Spontaneous binding of Hsc20 to Ssq1 and binding free energy

To study spontaneous binding between Hsc20 and Ssq1, the simulation dodecahedron box of size 16.2 nm × 16.2 nm × 16.2 nm containing single copies of Hsc20 and ATP-bound Ssq1 solvated with approx. 96500 water molecules was prepared. The proteins were initially placed such that the COM distance between positively charged residues of the Hsc20 J-domain helix II (R37, R38, R41) and negatively charged residues of Ssq1 NBD (D246, E248, D249, E253, D362, D364) was equal to 3 nm. For this system, five independent unbiased MD simulations were performed which were terminated after 60 ns as this time was sufficient for the native complex to spontaneously form in 2 out of 5 replicas. The free energy profiles for Ssq1 interaction with Hsc20 and with Hsc20(R37A,R41A) substitution variant were computed using replica-exchange umbrella sampling (REUS) [66]. As a reaction coordinate for REUS simulations, we used the COM distance between the respective binding sites corresponding to residues 34–44 of Hsc20 (Helix II of the J-domain) and residues 244–255 and 425–429 of Ssq1 (NBD). The harmonic potential was used to restrain the system along the reaction coordinate (spring constants and centers of the potential are summarized in S4 Table) and the initial configurations were taken from the spontaneous binding trajectories. The same initial configurations were used for the Hsc20(R37A,R41A) variant obtained by substituting R37 and R41 to alanine in all initial frames. For each of the REUS windows 400 ns long simulations were performed. Exchanges between neighboring windows were attempted every 2 ps and accepted by the Metropolis criterion. The free energy profiles were determined from the last 350 ns of thus obtained trajectories using the standard weighted histogram analysis method (WHAM) [63] and the Monte Carlo bootstrap method to estimate uncertainties of the free energy differences.

Evolutionary analyses

Protein sequences were retrieved from OMA orthology database [67]: Hsc20 (OMAGroup:684646), mitochondrial Hsp70 SSC1(OMAGroup:546796), DnaJ (OMAGroup:571688), DnaK (OMAGroup:546796) cytosolic Hsp70s SSA1(OMAGroup:555679), SSA2(OMAGroup:557653), SSA3(OMAGroup:556877), SSA4(OMAGroup:557493), endoplasmic reticulum Hsp70s KAR2 (OMAGroup:555433). Sequences were aligned using Clustal Omega v1.2.2 with default parameters [68]. Alignment logos were generated using the WebLogo server [69]. Ancestral states of sites homologous to D517 of Ssq1 and D481 of DnaK were reconstructed in FastML using empirical Bayes method with Maximum Likelihood reconstruction of insertions and deletions [70].

To infer Hsp70 phylogeny, a multiple sequence alignment was converted into a Hidden Markov Model [71] using hmmbuild program from the HMMER package. Forward–backward algorithm was used to compute a posterior probability (pp) for each site representing the degree of confidence in each position (residue or gap) of the alignment for each sequence. Amino acid positions with pp <0.7 were removed from the multiple sequence alignment [72]. 1,000 maximum likelihood (ML) searches were performed using RAxML v8.2.10 [73] with 1000 rapid bootstrap replicates, under the LG model of amino acid substitution and GAMMA model of rate heterogeneity with four discrete rate categories and the estimate of proportion of invariable sites (LG + I + G) [74], which was determined as the best-fit model by ProtTest v3.2 following Akaike Information Criterion [75]. Hsp70s tree topology was used for all calculations.

Maximum likelihood implementation of the Coev model [32] was used to test coevolution between amino acid sequence positions homologous to those occupied by residues that were at a cutoff distance of 0.5 nm across the Hsc20/Ssq1 and DnaJJD/DnaK binding interfaces. Coevolution analyses were based on the Hsp70 phylogeny, with Hsc20 orthologs paired with mitochondrial Hsp70 partners and DnaJ orthologs paired with DnaK partners. The fit of the Coev model in comparison to the 'independent' model was tested using Δ Akaike information criterion (ΔAIC = AICindependent—AICCoev). To establish the threshold for the ΔAIC value to be accepted as evidence for coevolution [32], we simulated sequence alignment based on the same phylogenetic tree as the original data but using the 'independent' model of sequence evolution. We used evolver software from the PAML 4 [76] for this simulation. The 99thpercentile of this expected ΔAIC distribution provided thresholds to consider the observed ΔAIC for Coev model to be accepted as coevolving with p<0.01 confidence level.

Protein purification

Expression of Bpa-containing Hsc20 variants (Q46Bpa and M51Bpa) with a polyhistidine tag at the C-terminus was performed using the E. coli strain C41(DE3) carrying two plasmids: pET21d encoding Hsc20 Bpa variant and pSUPT BpF encoding an engineered tRNA synthetase and tRNACUA for p-benzoyl-phenylalanine (Bpa) incorporation. Cells were cultured in the M9 minimal medium. Protein expression was induced by adding isopropyl-1-thio-D-galactopyranoside (IPTG) at a final concentration of 0.2 mM, L-arabinose at a final concentration of 0.2% (v/v) and Bpa at a final concentration of 0.1 mM. Cells were harvested and lysed using a French press in buffer J1 (20 mM Tris-HCl, pH 8.0, 500 mM NaCl, 1 mM phenylmethanesulfonylfluoride (PMSF), 10% glycerol, 30 mM imidazole, pH 8.0). After a clarifying spin, the supernatant was precipitated with ammonium sulfate (0.35 mg/ml). After centrifugation, the pellet was resuspended in buffer J1 and dialyzed overnight against buffer J1. Next, the proteins were subjected to HisBind Resin (Novagen) chromatography. After sequential washing steps with buffers J1, J2 (20 mM Tris-HCl, pH 8.0, 1 M NaCl, 1 mM PMSF, 10% glycerol, 30 mM imidazole, pH 8.0, 1 mM ATP, 2 mM MgCl2) and J1, proteins were eluted with a linear 30–300 mM imidazole gradient in buffer J1. Fractions containing Hsc20 were collected, pooled and dialyzed against buffer J3 (40 mM potassium phosphate, pH 6.8; 10% (v/v) glycerol; 75 mM NaCl; 0.05% Triton X-100; 5 mM β-mercaptoethanol). After dialysis proteins were loaded on a P11 cellulose column (Spark Scientific Ltd), washed with a 150 mM NaCl solution in J3 buffer and eluted with a linear 150–800 mM NaCl gradient in J3 buffer. Fractions containing the purest Hsc20 were collected, pooled and dialyzed against the final buffer (20 mM Tris-HCl pH 8.0, 10% glycerol, 50 mM NaCl; 5 mM β-mercaptoethanol). Aliquots were stored at –70°C.

All other Hsc20 variants with a polyhistidine tag at the C-terminus were purified according to [24], except E. coli strain C41(DE3) was used for expression. Ssq1 variants were purified as described [77]. Mge1 was purified as described [24]. Isu1 from Chaetomium thermophilum (Isu1Ct) was purified as described [25]

In all cases, protein concentrations, determined by using the Bradford (Bio-Rad) assay system with bovine serum albumin as a standard, are expressed as the concentration of monomers.

ATPase activity of Ssq1

The ATPase activity of Ssq1 variants defective in the Hsc20/Ssq1 interaction was measured as described by [24] with 0.5 μM Ssq1, 500 μM PVK-p peptide (LSLPPVKLHC) derived from the Isu1 sequence containing the PVK motif that interacts with substrate binding site of Ssq1 [78], 0.5 μM Mge1, and Hsc20 at the indicated concentrations in buffer A (40 mM Hepes–KOH, pH 7.5, 100 mM KCl, 1 mM dithiothreitol, 10 mM MgCl2, and 10% (v/v) glycerol). Reactions (15 μl) were initiated by the addition of ATP (2 μCi, DuPont NEG-003H, 3000 Ci/mmol) to a final concentration of 120 μM. Incubation was carried out at 25°C, and the reaction was terminated after 15 min by the addition of 100 μl of 1M perchloric acid and 1mM sodium phosphate. After addition of 400 μl of 20 mM ammonium molybdate and 400 μl of isopropyl acetate, samples were mixed and the phases were separated by a short centrifugation. An aliquot of the organic phase (150 μl), containing the radioactive orthophosphate-molybdate complex, was removed and radioactivity was determined by liquid scintillation counting. Control reactions lacking protein were included in all experiments. Values were plotted in GraphPadPrism. The ATPase activity corresponding to the maximal stimulation (Vm) of the wild-type Hsc20 was set to 1.

The ATPase activity of the Ssq1(R207D) variant was measured using an enzyme-coupled spectrophotometric assay [79]. Reaction mixtures of 0.5 ml contained 1 μM Ssq1 R207D or WT, 1 mM ATP, 0.265 mM NADH, 100 U/ml lactic dehydrogenase, 70 U/ml pyruvate kinase and 2.8 mM phosphoenolpyruvate in buffer (50 mM Hepes-KOH, pH 7.5, 150 mM KCl, 20 mM magnesium acetate, 10 mM dithiothreitol). If indicated, Mge1 was at 1 μM, Isu1Ct at 3 μM and Hsc20 at 3 μM. Reactions were initiated by mixing ATP with proteins and then incubating at 23°C. NADH absorbance was measured at 340 nm for 500 seconds in JASCO V-660 spectrophotometer.

Site specific photo-crosslinking

Reaction mixtures (20 μl) containing 20 μM Hsc20 Bpa-containing variants and 5 μM Ssq1 variants were prepared in the reaction buffer (40 mM Hepes–KOH, pH 7.5, 100 mM KCl, 1 mM dithiothreitol, 10 mM MgCl2, and 10% (v/v) glycerol). Reactions were initiated by the addition of ATP to a final concentration of 2 mM and incubated at 25°C for 15 min. Reactions were then exposed to 254 nm UV light (CL-1000 Ultraviolet Crosslinker) for 7 minutes. Control reactions were not exposed to UV light. Next, 5 μl of 4-fold concentrated Laemmli sample buffer (125 mM Tris-HCl, pH 6.8, 5% sodium dodecyl sulfate, 10% 2-mercaptoethanol, 20% (v/v) glycerol) was added to the reaction mixes. After 10 min incubation at 100°C, the total reaction volume was loaded onto SDS-PAGE and visualized by Coomassie blue stain. Every experiment was replicated three times.

Identification of photo-crosslinking products using mass spectroscopy (MS)

Gel pieces were dried with acetonitrile and subjected to reduction with 10 mM dithiothreitol in 100mM NH4HCO3 for 30 min at 57°C. Cysteines were then alkylated with 0.5 M iodoacetamide in 100 mM NH4HCO3 (45 minutes in a dark room at room temperature) and proteins were digested overnight with 10 ng/μl trypsin (CAT NO V5280, Promega) in 25 mM NH4HCO3 [80]. Digestion was stopped by adding trifluoroacetic acid at a final concentration of 0.1%. The mixture was centrifuged at 4°C, 14 000 g for 30 min, to remove any remaining solid. MS analysis was performed by LC-MS in the Laboratory of Mass Spectrometry (Institute of Biochemistry and Biophysics, Polish Academy of Sciences, Warsaw) using a nanoAcquity UPLC system (Waters) coupled to an Orbitrap Elite (Thermo Fisher Scientific). The mass spectrometer was operated in the data-dependent MS2 mode, and data were acquired in the m/z range of 300–2000. Peptides were separated by a 180 min linear gradient of 95% solution A (0.1% formic acid in water) to 35% solution B (acetonitrile and 0.1% formic acid). The measurement of each sample was preceded by three washing runs to avoid cross-contamination. Data were searched with Mascot [81], with the following parameters: enzyme: Trypsin, parent ions mass tolerance: 20 ppm, fragment ion mass tolerance: 0.1, missed cleavages: 1, fixed modifications: Carbamidomethyl (C), variable modifications: Oxidation (M). Protein identification was validated using a target-decoy search strategy.

Supporting information

S1 Fig.

Representative structures of Ssq1 (a) and Hsc20 (b) from our unbiased molecular dynamics simulations carried out using the CHARMM36 (left) and AMBER99SB-ILDN (right) force fields. Ssq1 and Hsc20 were simulated for 2 and 10 μs, respectively, with CHARMM, and for 2 and 5 μs, respectively, with AMBER. The obtained trajectories were subject to cluster analysis with an RMSD cutoff of 0.7 nm (Ssq1) and 0.35 nm (Hsc20). The centroids of each cluster were superimposed and are shown in different colors. NBD, SBDa and SBDb domains (in Ssq1) and J-domain and Isu1 binding domain (in Hsc20) are indicated.

https://doi.org/10.1371/journal.pcbi.1007913.s001

(PDF)

S2 Fig.

(a) Five most populated binding poses of Hsc20 along with their percentage contributions to the bound- state ensemble of the Hsc20-Ssq1 complex, obtained by the combined docking/MD simulation approach. Helices II and III are shown in cyan and the rest of Hsc20 in transparent cyan. SBD is in light brown, Ia and IIa subdomains of NBD are in light and dark grey, respectively, linker in yellow. (b) Time evolution of the average RMSD of helices II and III of the J-domain with respect to the center of the most populated binding pose above (the dominant bound state). Averaging was done over all 18 independent MD trajectories that were found to visit the dominant bound state at least once over the course of 500 ns. The shaded area shows the standard deviation of the RMSD.

https://doi.org/10.1371/journal.pcbi.1007913.s002

(PDF)

S3 Fig.

Contact maps of Ssq1 R207 and DnaK R167 calculated from 2D metadynamics trajectories in absence (grey bars) and presence (cyan bars) of Hsc20 and DnaJ J-domain, respectively. The contact is defined with a cutoff distance of 0.5 nm between non-hydrogen atoms, and the contacts with probability above 20% in either state are shown.

https://doi.org/10.1371/journal.pcbi.1007913.s003

(PDF)

S4 Fig. ATPase activity of Ssq1(R207D) variant.

Steady-state ATPase activity of Ssq1 was measured alone, or in the presence of indicated proteins, using an enzymatically coupled assay. Bars represent average value for three independent measurements with error bars as standard deviation. Concentration of components: 1 μM Ssq1; 3 μM Hsc20; 3 μM Isu1; 1 μM Mge1; 1 mM ATP.

https://doi.org/10.1371/journal.pcbi.1007913.s004

(PDF)

S5 Fig. Free energy profiles as a function of the distance d (defined in Fig 3) for Ssq1 and DnaK in red, for Hsc20-Ssq1 and DnaJJD-DnaK complexes in blue, and for Ssq1 and DnaK RNBD->D substitution variants in purple.

https://doi.org/10.1371/journal.pcbi.1007913.s005

(PDF)

S6 Fig. Phylogenetic distribution and ancestral reconstruction of aspartate and asparagine residues at the site homologous to DSBD of DnaK and Ssq1; amino acid state is indicated at the nodes of the maximum likelihood phylogeny of DnaK and mitochondrial Hsp70 orthologs from 289 species belonging to the indicated clades.

Scale is in amino acid substitutions per position.

https://doi.org/10.1371/journal.pcbi.1007913.s006

(PDF)

S7 Fig.

(a) Probabilities (p) of the most stable ion pairs across the J-domain/Ssq1 interface, calculated from the 10.5 μs trajectory of the dominant bound state. (b) Distance distributions of the ion pairs shown in (a).

https://doi.org/10.1371/journal.pcbi.1007913.s007

(PDF)

S8 Fig.

(a) Probabilities (p) of the most stable hydrogen bonds across the J-domain/Ssq1 interface, calculated from the 10.5 μs trajectory of the dominant bound state. (b) Relative energy inputs to Hsc20-Ssq1 complex binding energy calculated with MM/PBSA for the residues of the Hsc20 J-domain (top) and Ssq1 (bottom).

https://doi.org/10.1371/journal.pcbi.1007913.s008

(PDF)

S9 Fig. Identified products of Hsc20Bpa and Ssq1 WT site-specific UV-crosslinking.

(left) Purified Hsc20 Q46Bpa and M51Bpa variants, were incubated with Ssq1 WT in the presence of ATP, irradiated with UV light and separated by SDS-PAGE. Indicated bands (1–8) were excised and subjected to trypsin digestion followed by LC-MS. Migrations of size standard, in kDa, are indicated. (right) Identification of proteins detected in individual SDS-PAGE bands by Mascot, along with statistical significance of Peptide Mass Fingerprint search (Mascot score) and protein abundance estimation (emPAI).

https://doi.org/10.1371/journal.pcbi.1007913.s009

(PDF)

S10 Fig. Site specific crosslinking of Hsc20 to Ssq1.

Purified Hsc20 S34Bpa, Q46Bpa and M51Bpa variants, as indicated, were incubated with purified Ssq1 in the presence of ATP. After UV irradiation (+), or as a control no irradiation (-), reaction mixtures were separated by SDS-PAGE. Migrations of size standard, in kDa, are indicated; Q:J* and Q:J** indicate positions of the Ssq1-Hsc20 crosslinks, which were identified using mass-spectroscopy (S9 Fig).

https://doi.org/10.1371/journal.pcbi.1007913.s010

(PDF)

S11 Fig.

Hsc20/Ssq1 recognition is driven by long-range electrostatic forces (a) Electrostatic isopotential contours at +/- 1kT/e (blue and red, respectively) around Hsc20 and Ssq1 and (b) around Hsc20 double substitution variant R37A,R41A. Coloring of structures as in Fig 1: Hsc20: J-domain (cyan); Isu binding domain (light blue); Ssq1: SBDα and SBDβ (brown); interdomain linker (yellow); NBD subdomains Ia and IIa in light and dark grey, respectively. (c) Structural alignment of Hsc20-Ssq1 complex obtained by spontaneous binding MD simulations and the dominant bound state obtained by molecular docking/MD simulations. The heavy-atom RMSD for J-domains after aligning NBDs equals to 2.5 Å. (d) Free energy profiles for Hsc20-Ssq1 binding as a function of separation distance, d, between J-domain of Hsc20 and NBD subdomain IIa of Ssq1. Binding of Hsc20 wild-type (WT) (blue line), binding of Hsc20(R37A,R41A) variant (red line); shaded areas show the standard error. The vertical dashed line, demarcating the bund and unbound states, indicates the distance beyond which all specific interactions between the binding partners are lost. Representative structures are shown in circles: unbound state (left), bound state—corresponding to the energy well (right).

https://doi.org/10.1371/journal.pcbi.1007913.s011

(PDF)

S12 Fig.

Sequence variability of positions that form the Hsc20/Ssq1 binding interfaces (positions shown in green in Fig 5A) for orthologs of Hsc20 (left) and orthologs of Ssq1 (right) from the following taxonomic groups: Saccharomycotina (true yeasts), Pezizomycotina (p. yeasts), Taphiromycotina (t. yeasts), Basidiomycota (club fungi), Animals (animals) and Plants (plants). Sequence logos represent the amino acid frequency of each interfacial position in the orthologs examined, with positively charged (blue), negatively charged (red) and uncharged (grey) residues. Position numbering is that for S. cerevisiae. Phylogenetic relationships among taxonomic groups are depicted on the left.

https://doi.org/10.1371/journal.pcbi.1007913.s012

(PDF)

S13 Fig.

Sequence variability of positions that form the DnaJ/DnaK binding interfaces (positions shown in green in Fig 5A) for orthologs of DnaJ (top) and orthologs of DnaK (bottom) from bacteria and from mitochondria for the following taxonomic groups: Bacteria, Saccharomycotina (true yeasts), Pezizomycotina (p. yeasts), Taphiromycotina (t. yeasts), Basidiomycota (club fungi), Animals (animals) and Plants (plants). Sequence logos represent the amino acid frequency of each interfacial position in the orthologs examined, with positively charged (blue), negatively charged (red) and uncharged (grey) residues. Position numbering is that for E. coli. Phylogenetic relationships among taxonomic groups are depicted on the left.

https://doi.org/10.1371/journal.pcbi.1007913.s013

(PDF)

S14 Fig. Maximum likelihood phylogeny of mitochondrial orthologs of Ssq1 from 67 Ascomycota species used in the sequence variability and co-evolution analyses.

Major clades are indicated.

https://doi.org/10.1371/journal.pcbi.1007913.s014

(PDF)

S15 Fig. Maximum likelihood phylogeny of 117 bacterial DnaK orthologs used in the sequence variability and coevolution analyses.

Major clades are indicated.

https://doi.org/10.1371/journal.pcbi.1007913.s015

(PDF)

S16 Fig. ΔAIC (AICindependent—AICCoev) distribution calculated for each pair of positions derived from amino acid sequence of length 33 simulated under the independent evolution model using mitochondrial Hsp70 phylogeny (S12 Fig).

The co-evolution acceptance threshold with p>0.01 confidence was established as ΔAIC = 7.38 (red dashed line).

https://doi.org/10.1371/journal.pcbi.1007913.s016

(PDF)

S17 Fig. ΔAIC (AICindependent—AICCoev) distribution calculated for each pair of positions derived from amino acid sequence of length 54 simulated under the independent evolution model using DnaK phylogeny (S13 Fig).

The co-evolution acceptance threshold with p>0.01 confidence was established as ΔAIC = 9.473 (red dashed line).

https://doi.org/10.1371/journal.pcbi.1007913.s017

(PDF)

S1 Table. Probability of observing RNBD in four different states.

https://doi.org/10.1371/journal.pcbi.1007913.s018

(PDF)

S2 Table. Statistical support for pairs of coevolving Hsc20/Hsp70 positions from Ascomycota.

https://doi.org/10.1371/journal.pcbi.1007913.s019

(PDF)

S3 Table. Statistical support for pairs of coevolving DnaJJD/DnaK positions from bacteria.

https://doi.org/10.1371/journal.pcbi.1007913.s020

(PDF)

S4 Table. Force constants used in Umbrella Sampling for Hsc20-Ssq1 spontaneous binding.

https://doi.org/10.1371/journal.pcbi.1007913.s021

(PDF)

S1 Movie. Spontaneous transition of RNBD to the intermolecular contact with DHPD.

https://doi.org/10.1371/journal.pcbi.1007913.s022

(MP4)

S2 Movie. Spontaneous Hsc20-Ssq1 complex formation in MD.

https://doi.org/10.1371/journal.pcbi.1007913.s023

(MP4)

Acknowledgments

We thank Xavier Meyer and Daniele Silvestro, University of Lausanne, for discussions about detecting coevolution. We thank Aneta Grabinska for technical help with the ATPase analyses.

References

  1. 1. Kampinga HH, Craig EA. The HSP70 chaperone machinery: J proteins as drivers of functional specificity. Nature Reviews Molecular Cell Biology. 2010;11(8):579–92. pmid:20651708; PubMed Central PMCID: PMC3003299.
  2. 2. Rosenzweig R, Nillegoda NB, Mayer MP, Bukau B. The Hsp70 chaperone network. Nature Reviews Molecular Cell Biology. 2019;353:aac4354. pmid:31253954.
  3. 3. Clerico EM, Meng W, Pozhidaeva A, Bhasne K, Petridis C, Gierasch LM. Hsp70 molecular chaperones: multifunctional allosteric holding and unfolding machines. Biochemical Journal. 2019;476(11):1653–77. pmid:31201219.
  4. 4. Craig EA, Marszalek J. How Do J-Proteins Get Hsp70 to Do So Many Different Things? Trends in Biochemical Sciences. 2017;42(5):355–68. pmid:28314505; PubMed Central PMCID: PMC5409888.
  5. 5. Kityk R, Kopp J, Sinning I, Mayer MP. Structure and Dynamics of the ATP-Bound Open Conformation of Hsp70 Chaperones. Molecular Cell. 2012;48(6):1–12. pmid:23123194
  6. 6. Qi R, Sarbeng EB, Liu Q, Le KQ, Xu X, Xu H, et al. Allosteric opening of the polypeptide-binding site when an Hsp70 binds ATP. Nature Structural & Molecular Biology. 2013;20(7):1–10. pmid:23708608
  7. 7. Mapa K, Sikor M, Kudryavtsev V, Waegemann K, Kalinin S, Seidel CAM, et al. The conformational dynamics of the mitochondrial Hsp70 chaperone. Molecular Cell. 2010;38(1):89–100. pmid:20385092.
  8. 8. Marcinowski M, Höller M, Feige MJ, Baerend D, Lamb DC, Buchner J. Substrate discrimination of the chaperone BiP by autonomous and cochaperone-regulated conformational transitions. Nature Structural & Molecular Biology. 2011;18(2):150–8. pmid:21217698.
  9. 9. Banerjee R, Jayaraj GG, Peter JJ, Kumar V, Mapa K. Monitoring conformational heterogeneity of the lid of DnaK substrate-binding domain during its chaperone cycle. The FEBS journal. 2016;283(15):2853–68. pmid:27248857.
  10. 10. Zhuravleva A, Clerico EM, Gierasch LM. An interdomain energetic tug-of-war creates the allosterically active state in Hsp70 molecular chaperones. Cell. 2012;151(6):1296–307. pmid:23217711; PubMed Central PMCID: PMC3521165.
  11. 11. Kityk R, Vogel M, Schlecht R, Bukau B, Mayer MP. Pathways of allosteric regulation in Hsp70 chaperones. Nature Communications. 2015;6:1–11. 16370000797795331459related:g61zXf7sLeMJ. pmid:26383706
  12. 12. Vogel M, Bukau B, Mayer MP. Allosteric Regulation of Hsp70 Chaperones by a Proline Switch. Molecular Cell. 2006;21(3):359–67. pmid:16455491
  13. 13. General IJ, Liu Y, Blackburn ME, Mao W, Gierasch LM, Bahar I. ATPase Subdomain IA Is a Mediator of Interdomain Allostery in Hsp70 Molecular Chaperones. PLoS Computational Biology. 2014;10(5):e1003624. pmid:24831085
  14. 14. Vogel M, Mayer MP, Bukau B. Allosteric regulation of Hsp70 chaperones involves a conserved interdomain linker. Journal of Biological Chemistry. 2006;281(50):38705–11. pmid:17052976.
  15. 15. Swain JF, Dinler G, Sivendran R, Montgomery DL, Stotz M, Gierasch LM. Hsp70 Chaperone Ligands Control Domain Association via an Allosteric Mechanism Mediated by the Interdomain Linker. Molecular Cell. 2007;26(1):27–39. pmid:17434124
  16. 16. Alderson TR, Kim JH, Cai K, Frederick RO, Tonelli M, Markley JL. The Specialized Hsp70 (HscA) Interdomain Linker Binds to Its Nucleotide-Binding Domain and Stimulates ATP Hydrolysis in Both cis and trans Configurations. Biochemistry. 2014;53(46):7148–59. pmid:25372495
  17. 17. English CA, Sherman W, Meng W, Gierasch LM. The Hsp70 interdomain linker is a dynamic switch that enables allosteric communication between two structured domains. Journal of Biological Chemistry. 2017;292(36):14765–74. pmid:28754691; PubMed Central PMCID: PMC5592658.
  18. 18. Lai AL, Clerico EM, Blackburn ME, Patel NA, Robinson CV, Borbat PP, et al. Key features of an Hsp70 chaperone allosteric landscape revealed by ion-mobility native mass spectrometry and double electron-electron resonance. Journal of Biological Chemistry. 2017;292(21):8773–85. pmid:28428246.
  19. 19. Liberek K, Marszalek J, Ang D, Georgopoulos C, Zylicz M. Escherichia coli DnaJ and GrpE heat shock proteins jointly stimulate ATPase activity of DnaK. Proceedings of the National Academy of Sciences of the United States of America. 1991:1–5.
  20. 20. Suh WC, Burkholder WF, Lu CZ, Zhao X, Gottesman ME, Gross CA. Interaction of the Hsp70 molecular chaperone, DnaK, with its cochaperone DnaJ. Proceedings of the National Academy of Sciences of the United States of America. 1998;95(26):15223. 10894679841592021999related:75Oz65etMZcJ. pmid:9860950
  21. 21. Gassler CS, Buchberger A, Laufen T, Mayer MP, Schröder H, Valencia A, et al. Mutations in the DnaK chaperone affecting interaction with the DnaJ cochaperone. Proceedings of the National Academy of Sciences of the United States of America. 1998;95(26):15229–34. pmid:9860951; PubMed Central PMCID: PMC28025.
  22. 22. Kityk R, Kopp J, Mayer MP. Molecular Mechanism of J-Domain-Triggered ATP Hydrolysis by Hsp70 Chaperones. Molecular Cell. 2018;69(2):227–37.e4. pmid:29290615.
  23. 23. Malinverni D, Jost Lopez A, De Los Rios P, Hummer G, Barducci A. Modeling Hsp70/Hsp40 interaction by multi-scale molecular simulations and coevolutionary sequence analysis. eLife. 2017;6:19. pmid:28498104; PubMed Central PMCID: PMC5519331.
  24. 24. Dutkiewicz R, Schilke B, Knieszner H, Walter W, Craig EA, Marszalek J. Ssq1, a mitochondrial Hsp70 involved in iron-sulfur (Fe/S) center biogenesis. Similarities to and differences from its bacterial counterpart. Journal of Biological Chemistry. 2003;278(32):29719–27. pmid:12756240.
  25. 25. Delewski W, Paterkiewicz B, Manicki M, Schilke B, Tomiczek B, Ciesielski SJ, et al. Iron–Sulfur Cluster Biogenesis Chaperones: Evidence for Emergence of Mutational Robustness of a Highly Specific Protein–Protein Interaction. Molecular Biology and Evolution. 2016;33(3):643–56. pmid:26545917
  26. 26. Dror RO, Dirks RM, Grossman JP, Xu H, Shaw DE. Biomolecular Simulation: A Computational Microscope for Molecular Biology. Annual Review of Biophysics. 2012;41(1):429–52. pmid:22577825
  27. 27. Barducci A, Bonomi M, Parrinello M. Metadynamics. Wiley Online Library. 2011. 264BBD0C-65CE-4ABD-B971-C2D692E79774.
  28. 28. Valsson O, Tiwary P, Parrinello M. Enhancing Important Fluctuations: Rare Events and Metadynamics from a Conceptual Viewpoint. Annual Review of Physical Chemistry. 2016;67(1):159–84. pmid:26980304
  29. 29. Donald JE, Kulp DW, DeGrado WF. Salt bridges: geometrically specific, designable interactions. Proteins: Structure, Function, and Bioinformatics. 2011;79(3):898–915. pmid:21287621; PubMed Central PMCID: PMC3069487.
  30. 30. Krishnamurthy M, Dugan A, Nwokoye A, Fung Y-H, Lancia JK, Majmudar CY, et al. Caught in the act: covalent cross-linking captures activator-coactivator interactions in vivo. ACS Chemical Biology. 2011;6(12):1321–6. pmid:21977905; PubMed Central PMCID: PMC3245988.
  31. 31. de Juan D, Pazos F, Valencia A. Emerging methods in protein co-evolution. Nature Reviews Genetics. 2013;14(4):1–13. pmid:23458856
  32. 32. Dib L, Silvestro D, Salamin N. Evolutionary footprint of coevolving positions in genes. Bioinformatics. 2014;30(9):1241–9. Epub 2nd edn. pmid:24413673.
  33. 33. Rist W, Graf C, Bukau B, Mayer MP. Amide hydrogen exchange reveals conformational changes in hsp70 chaperones important for allosteric regulation. Journal of Biological Chemistry. 2006;281(24):16493–501. pmid:16613854.
  34. 34. Nicolaï A, Delarue P, Senet P. Decipher the mechanisms of protein conformational changes induced by nucleotide binding through free-energy landscape analysis: ATP binding to Hsp70. PLoS Computational Biology. 2013;9(12):e1003379. pmid:24348227; PubMed Central PMCID: PMC3861046.
  35. 35. Jiang J, Maes EG, Taylor AB, Wang L, Hinck AP, Lafer EM, et al. Structural Basis of J Cochaperone Binding and Regulation of Hsp70. Molecular Cell. 2007;28(3):422–33. pmid:17996706
  36. 36. Zhuravleva A, Gierasch LM. Substrate-binding domain conformational dynamics mediate Hsp70 allostery. Proceedings of the National Academy of Sciences of the United States of America. 2015;112(22):E2865–E73. pmid:26038563
  37. 37. Sikor M, Mapa K, von Voithenberg LV, Mokranjac D, Lamb DC. Real-time observation of the conformational dynamics of mitochondrial Hsp70 by spFRET. The EMBO Journal. 2013:1–11.
  38. 38. Wieteska L, Shahidi S, Zhuravleva A. Allosteric fine-tuning of the conformational equilibrium poises the chaperone BiP for post-translational regulation. eLife. 2017;6:18966. pmid:29064369.
  39. 39. Meng W, Clerico EM, McArthur N, Gierasch LM. Allosteric landscapes of eukaryotic cytoplasmic Hsp70s are shaped by evolutionary tuning of key interfaces. Proceedings of the National Academy of Sciences of the United States of America. 2018;115(47):11970–5. pmid:30397123; PubMed Central PMCID: PMC6255180.
  40. 40. Stetz G, Verkhivker GM. Computational Analysis of Residue Interaction Networks and Coevolutionary Relationships in the Hsp70 Chaperones: A Community-Hopping Model of Allosteric Regulation and Communication. PLoS Computational Biology. 2017;13(1):e1005299. pmid:28095400; PubMed Central PMCID: PMC5240922.
  41. 41. Bascos NAD, Mayer MP, Bukau B, Landry SJ. The Hsp40 J-domain modulates Hsp70 conformation and ATPase activity with a semi-elliptical spring. Protein Science. 2017;26(9):1838–51. pmid:28685898
  42. 42. Schilke B, Williams B, Knieszner H, Pukszta S, D&apos;Silva P, Craig EA, et al. Evolution of Mitochondrial Chaperones Utilized in Fe-S Cluster Biogenesis. Current Biology. 2006;16(16):1660–5. pmid:16920629
  43. 43. Pukszta S, Schilke B, Dutkiewicz R, Kominek J, Moczulska K, Stepien B, et al. Co-evolution-driven switch of J-protein specificity towards an Hsp70 partner. EMBO Reports. 2010;11(5):360–5. pmid:20224575
  44. 44. Meyer AE, Hung N-J, Yang P, Johnson AW, Craig EA. The specialized cytosolic J-protein, Jjj1, functions in 60S ribosomal subunit biogenesis. Proceedings of the National Academy of Sciences of the United States of America. 2007;104(5):1558–63. pmid:17242366; PubMed Central PMCID: PMC1785244.
  45. 45. Sahi C, Craig EA. Network of general and specialty J protein chaperones of the yeast cytosol. Proceedings of the National Academy of Sciences of the United States of America. 2007;104(17):7163–8. pmid:17438278; PubMed Central PMCID: PMC1855418.
  46. 46. Ahmad A, Bhattacharya A, McDonald RA, Cordes M, Ellington B, Bertelsen EB, et al. Heat shock protein 70 kDa chaperone/DnaJ cochaperone complex employs an unusual dynamic interface. Proceedings of the National Academy of Sciences of the United States of America. 2011;108(47):18966–71. pmid:22065753; PubMed Central PMCID: PMC3223468.
  47. 47. Sousa R, Jiang J, Lafer EM, Lafer EM, Hinck AP, Hinck AP, et al. Evaluation of competing J domain:Hsp70 complex models in light of existing mutational and NMR data. Proceedings of the National Academy of Sciences of the United States of America. 2012;109(13):E734–author reply E5. pmid:22403069; PubMed Central PMCID: PMC3323956.
  48. 48. Garimella R, Liu X, Qiao W, Liang X, Zuiderweg ERP, Riley MI, et al. Hsc70 Contacts Helix III of the J Domain from Polyomavirus T Antigens: Addressing a Dilemma in the Chaperone Hypothesis of How They Release E2F from pRb †. Biochemistry. 2006;45(22):6917–29. pmid:16734427
  49. 49. Zarouchlioti C, Parfitt DA, Li W, Gittings LM, Cheetham ME. DNAJ Proteins in neurodegeneration: essential and protective factors. Philosophical transactions of the Royal Society of London Series B, Biological Sciences. 2017;373(1738):20160534–14. pmid:29203718
  50. 50. Calderwood SK, Gong J. Heat Shock Proteins Promote Cancer: It&apos;s a Protection Racket. Trends in Biochemical Sciences. 2016;41(4):311–23. pmid:26874923
  51. 51. Li X, Shao H, Taylor IR, Gestwicki JE. Targeting Allosteric Control Mechanisms in Heat Shock Protein 70 (Hsp70). Current Topics in Medicinal Chemistry. 2016;16(25):2729–40. pmid:27072701; PubMed Central PMCID: PMC5502483.
  52. 52. Abraham MJ, Murtola T, Schulz R, Páll S, Smith JC, Hess B, et al. GROMACS: High performance molecular simulations through multi-level parallelism from laptops to supercomputers. SoftwareX. 2015;1–2:19–25.
  53. 53. Bonomi M, Branduardi D, Bussi G, Camilloni C, Provasi D, Raiteri P, et al. PLUMED: a portable plugin for free-energy calculations with molecular dynamics 2009:1–15.
  54. 54. Huang J, MacKerell AD Jr. CHARMM36 all-atom additive protein force field: Validation based on comparison to NMR data. Journal of Computational Chemistry. 2013;34(25):2135–45. pmid:23832629
  55. 55. Bussi G, Donadio D, Parrinello M. Canonical sampling through velocity rescaling. The Journal of Chemical Physics. 2007;126(1):014101–8. Epub 3. pmid:17212484
  56. 56. Parrinello M, Rahman A. Polymorphic transitions in single crystals: A new molecular dynamics method. Journal of Applied Physics. 1981;52(12):7182–90.
  57. 57. Darden T, York D, Pedersen L. Particle mesh Ewald: An N⋅log (N) method for Ewald sums in large systems. Journal of Chemical Physics. 1993;98(12):10089–92.
  58. 58. Yang J, Zhang Y. Protein Structure and Function Prediction Using I-TASSER. Current Protocols in Bioinformatics. 2015;52:5.8.1–15. pmid:26678386; PubMed Central PMCID: PMC4871818.
  59. 59. Lindorff-Larsen K, Piana S, Palmo K, Maragakis P, Klepeis JL, Dror RO, et al. Improved side-chain torsion potentials for the Amber ff99SB protein force field. Proteins: Structure, Function, and Bioinformatics. 2010;19:NA-NA. pmid:20408171
  60. 60. Kozakov D, Hall DR, Xia B, Porter KA, Padhorny D, Yueh C, et al. The ClusPro web server for protein-protein docking. Nature Protocols. 2017;12(2):255–78. pmid:28079879; PubMed Central PMCID: PMC5540229.
  61. 61. Raiteri P, Laio A, Gervasio FL, Micheletti C, Parrinello M. Efficient Reconstruction of Complex Free Energy Landscapes by Multiple Walkers Metadynamics †. Journal of Physical Chemistry B. 2006;110(8):3533–9. pmid:16494409
  62. 62. Laio A, Gervasio FL. Metadynamics: a method to simulate rare events and reconstruct the free energy in biophysics, chemistry and material science. Reports on Progress in Physics. 2008;71(12):126601–23.
  63. 63. Kumar S, Bouzida D, Swendsen RH, Kollman PA, Rosenberg JM. The Weighted Histogram Analysis Method for Free-Energy Calculations on Biomolecules. Journal of Computational Chemistry. 1992;13:1011–21.
  64. 64. Kumar S, Nussinov R. Relationship between ion pair geometries and electrostatic strengths in proteins. Biophysical Journal. 2002;83(3):1595–612. pmid:12202384; PubMed Central PMCID: PMC1302257.
  65. 65. Kumari R, Kumar R, Consortium OSDD, Lynn A. g_mmpbsa—A GROMACS Tool for High-Throughput MM-PBSA Calculations. Journal of Chemical Information and Modeling. 2014;54(7):1951–62. pmid:24850022
  66. 66. Sugita Y, Kitao A, Okamoto Y. Multidimensional replica-exchange method for free-energy calculations. Journal of Chemical Physics. 2000;113(15):6042–51.
  67. 67. Altenhoff AM, Glover NM, Train C-M, Kaleb K, Warwick Vesztrocy A, Dylus D, et al. The OMA orthology database in 2018: retrieving evolutionary relationships among all domains of life through richer web and programmatic interfaces. Nucleic Acids Research. 2018;46(D1):D477–D85. pmid:29106550; PubMed Central PMCID: PMC5753216.
  68. 68. Sievers F, Wilm A, Dineen D, Gibson TJ, Karplus K, Li W, et al. Fast, scalable generation of high-quality protein multiple sequence alignments using Clustal Omega. Molecular Systems Biology. 2011;7:1–6. pmid:21988835
  69. 69. Crooks GE, Hon G, Chandonia J-M, Brenner SE. WebLogo: a sequence logo generator. Genome Research. 2004;14(6):1188–90. pmid:15173120; PubMed Central PMCID: PMC419797.
  70. 70. Ashkenazy H, Penn O, Doron-Faigenboim A, Cohen O, Cannarozzi G, Zomer O, et al. FastML: a web server for probabilistic reconstruction of ancestral sequences. Nucleic Acids Research. 2012;40(W1):W580–W4. pmid:22661579
  71. 71. Eddy SR. Profile hidden Markov models. Bioinformatics. 1998;14(9):755–63. pmid:9918945
  72. 72. Sarangi GK, Romagné F, Castellano S. Distinct Patterns of Selection in Selenium-Dependent Genes between Land and Aquatic Vertebrates. Molecular Biology and Evolution. 2018;35(7):1744–56. pmid:29669130
  73. 73. Stamatakis A. RAxML version 8: a tool for phylogenetic analysis and post-analysis of large phylogenies. Bioinformatics. 2014;30(9):1312–3. pmid:24451623; PubMed Central PMCID: PMC3998144.
  74. 74. Le SQ, Gascuel O. An Improved General Amino Acid Replacement Matrix. Molecular Biology and Evolution. 2008;25(7):1307–20. pmid:18367465
  75. 75. Darriba D, Taboada GL, Doallo R, Posada D. ProtTest 3: fast selection of best-fit models of protein evolution. Bioinformatics. 2011;27(8):1164–5. pmid:21335321
  76. 76. Yang Z. PAML 4: Phylogenetic Analysis by Maximum Likelihood. Molecular Biology and Evolution. 2007;24(8):1586–91. pmid:17483113
  77. 77. Manicki M, Majewska J, Ciesielski S, Schilke B, Blenska A, Kominek J, et al. Overlapping Binding Sites of the Frataxin Homologue Assembly Factor and the Heat Shock Protein 70 Transfer Factor on the Isu Iron-Sulfur Cluster Scaffold Protein. Journal of Biological Chemistry. 2014;289(44):30268–78. pmid:25228696
  78. 78. Knieszner H, Schilke B, Dutkiewicz R, D&apos;Silva P, Cheng S, Ohlson M, et al. Compensation for a defective interaction of the hsp70 ssq1 with the mitochondrial Fe-S cluster scaffold isu. Journal of Biological Chemistry. 2005;280(32):28966–72. pmid:15958384.
  79. 79. Nørby JG. Coupled assay of Na+,K+-ATPase activity. Methods in Enzymology. 1988;156:116–9. pmid:2835597.
  80. 80. Orlowska KP, Klosowska K, Szczesny RJ, Cysewski D, Krawczyk PS, Dziembowski A. A new strategy for gene targeting and functional proteomics using the DT40 cell line. Nucleic Acids Research. 2013;41(17):e167. pmid:23892402; PubMed Central PMCID: PMC3783193.
  81. 81. Perkins DN, Pappin DJ, Creasy DM, Cottrell JS. Probability-based protein identification by searching sequence databases using mass spectrometry data. Electrophoresis. 1999;20(18):3551–67. Epub 1999/12/28. pmid:10612281.