Probing Allosteric Inhibition Mechanisms of the ... - ACS Publications

Jul 22, 2016 - Technology, Chapman University, One University Drive, Orange, California 92866, United States. ‡. Chapman University School of Pharma...
1 downloads 0 Views 9MB Size
Subscriber access provided by CORNELL UNIVERSITY LIBRARY

Article

Probing Allosteric Inhibition Mechanisms of the Hsp70 Chaperone Proteins Using Molecular Dynamics Simulations and Analysis of the Residue Interaction Networks Gabrielle Stetz, and Gennady M. Verkhivker J. Chem. Inf. Model., Just Accepted Manuscript • DOI: 10.1021/acs.jcim.5b00755 • Publication Date (Web): 22 Jul 2016 Downloaded from http://pubs.acs.org on July 27, 2016

Just Accepted “Just Accepted” manuscripts have been peer-reviewed and accepted for publication. They are posted online prior to technical editing, formatting for publication and author proofing. The American Chemical Society provides “Just Accepted” as a free service to the research community to expedite the dissemination of scientific material as soon as possible after acceptance. “Just Accepted” manuscripts appear in full in PDF format accompanied by an HTML abstract. “Just Accepted” manuscripts have been fully peer reviewed, but should not be considered the official version of record. They are accessible to all readers and citable by the Digital Object Identifier (DOI®). “Just Accepted” is an optional service offered to authors. Therefore, the “Just Accepted” Web site may not include all articles that will be published in the journal. After a manuscript is technically edited and formatted, it will be removed from the “Just Accepted” Web site and published as an ASAP article. Note that technical editing may introduce minor changes to the manuscript text and/or graphics which could affect content, and all legal disclaimers and ethical guidelines that apply to the journal pertain. ACS cannot be held responsible for errors or consequences arising from the use of information contained in these “Just Accepted” manuscripts.

Journal of Chemical Information and Modeling is published by the American Chemical Society. 1155 Sixteenth Street N.W., Washington, DC 20036 Published by American Chemical Society. Copyright © American Chemical Society. However, no copyright claim is made to original U.S. Government works, or works produced by employees of any Commonwealth realm Crown government in the course of their duties.

Page 1 of 114

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Information and Modeling

Probing Allosteric Inhibition Mechanisms of the Hsp70 Chaperone Proteins Using Molecular Dynamics Simulations and Analysis of the Residue Interaction Networks Gabrielle Stetz1, Gennady M. Verkhivker1, 2

1



Graduate Program in Computational and Data Sciences, Department of Computational Sciences,

Schmid College of Science and Technology, Chapman University, One University Drive, Orange, CA 92866, USA 2

Chapman University School of Pharmacy, Irvine, CA 92618, USA

‡corresponding author E-mail: [email protected]

ACS Paragon Plus Environment

1

Journal of Chemical Information and Modeling

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 2 of 114

Abstract Although molecular mechanisms of allosteric regulation in the Hsp70 chaperones have been extensively studied at both structural and functional levels, allosteric study,

the current understanding of

inhibition of chaperone activities by small molecules is still lacking. In the current using a battery of computational approaches we

probed

allosteric inhibition

mechanisms of E. coli Hsp70 (DnaK) and human Hsp70 proteins by small molecule inhibitors PET-16 and novolactone. were combined with

Molecular dynamics simulations and binding free energy analysis network-based modeling

of

residue interactions and allosteric

communications to systematically characterize and compare molecular signatures of the apo form, substrate-bound, and inhibitor-bound chaperone complexes. mechanism by which the allosteric inhibitors may

The results suggested a

leverage binding energy hotspots in the

interaction networks to stabilize a specific conformational state and impair the inter-domain allosteric control.

Using the network-based centrality analysis and community detection, we

demonstrated that substrate binding may strengthen the connectivity of local interaction communities, leading to a dense interaction network that can promote communication.

an efficient allosteric

In contrast, binding of PET-16 to DnaK may induce significant dynamic

changes and lead to a fractured interaction network and impaired allosteric communications in the DnaK complex. By using a mechanistic-based analysis of distance fluctuation maps and allosteric propensities of protein residues, we determined that the allosteric network in the PET16 complex

may be

small and localized due to the reduced communication and low

cooperativity of the substrate binding loops, which may promote the higher rates of substrate dissociation and the decreased substrate affinity. In comparison with the significant PET-16, binding of novolactone to HSPA1A may cause

effect of

only moderate network changes and

ACS Paragon Plus Environment

2

Page 3 of 114

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Information and Modeling

preserve allosteric coupling between the allosteric pocket and the substrate binding region. The impact of novolactone on the conformational dynamics and allosteric communications in the HSPA1A complex was comparable to the substrate effect, which is consistent with the experimental evidence that PET-16, but not novolactone binding, can significantly decrease substrate affinity. We argue that the unique dynamic and network signatures of PET-16 and novolactone may be linked with the experimentally observed

functional effects

of these

inhibitors on allosteric regulation and substrate binding.

Introduction

The mammalian heat-inducible cytosolic 70-kilodalton (kDa) Heat shock protein Hsp70 along with its evolutionarily conserved bacterial orthologue, DnaK belong to a conserved and ubiquitous

chaperone family that is involved in regulation of diverse cellular processes,

including heat shock response and maintenance of protein stability.1-14 The functions of the Hsp70 chaperones in the protein folding networks and their remarkable adaptability to various substrates are regulated through an allosteric mechanism that involves cooperative involvement of the nucleotide, substrate, co-chaperones and nucleotide exchange factor proteins.15,16 functional cycle of the Hsp70 chaperones is

The

mainly driven by the nucleotide-dependent

allosteric transitions between ATP-bound and ADP-bound chaperone states. In the course of this cycle, binding of ATP in the nucleotide-binding domain (NBD) can allosterically induce the increased rates of substrate dissociation and lower the substrate binding affinity in the substratebinding domain (SBD).

Structural understanding of the Hsp70 mechanisms has been steadily

advancing in recent years,

evolving from the initial structural studies of the isolated domains

ACS Paragon Plus Environment

3

Journal of Chemical Information and Modeling

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 4 of 114

to high-resolution structures of the full-length two-domain constructs and complexes with the nucleotides and small molecule inhibitors. Crystal structures and nuclear magnetic resonance (NMR) studies of the NBD constructs17-20 and SBD complexes with the peptide substrates21-24 have provided first important insights into the structural basis of substrate specificity with Hsp70.

Amide hydrogen exchange and mass spectrometry (HX-MS) studies of the two-

domain

DnaK constructs25-28

along with

subsequent biophysical experiments29,30 have

shown that the NBD and SBD units are largely independent in a domain-undocked form of the ADP-bound and nucleotide-free states, whereas ATP binding could stabilize a domain-docked structure that facilitates the inter-domain signaling. The solution NMR structure of the fulllength DnaK has confirmed that the NBD and SBD are spatially separated and only weakly bound in the ADP-bound state.31 Crystal structures of the ATP-bound DnaK constructs32,33 have revealed the molecular details of the extensive NBD-SBD interface in the domaindocked

form.

Electron paramagnetic resonance spectroscopy34

and

single molecule

fluorescence experiments35,36 have confirmed the highly dynamic nature of both ADP-bound and ATP-bound DnaK structures, revealing significant movements of the SBD-α helical subdomain (often termed as a helical lid) and presence of multiple conformational states. Thermodynamics parameters of

functional DnaK states obtained

in isothermal titration

calorimetry (ITC) experiments37 have suggested a dynamic equilibrium between ATP-bound states and ADP-bound chaperone conformations that is accompanied by switching between low and high substrate affinity. NMR studies of allosteric regulation in DnaK have shown that

ATP binding could promote conformational flexibility of the substrate binding loops,

leading to the low affinity and fast kinetics of substrate release and association.38 This study

ACS Paragon Plus Environment

4

Page 5 of 114

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Information and Modeling

has proposed that structural changes and allosteric communication between the substrate binding loops and the SBD-NBD interface may be central to the allosteric mechanism in DnaK.

Site-directed mutagenesis studies have

been widely employed to examine role of specific

functional residues in allosteric regulation of DnaK.39-44 A number of functionally important sites are located in the SBD regions, and N415G, P419A and D326V)

may have

mutations of these residues ( for instance K414I, significant detrimental effect on

allostery.41,44

Functional analysis of DnaK-SBD mutants has confirmed that the major structural elements responsible

for the thermodynamics and kinetics of substrate binding include the central

hydrophobic pocket, the helical lid and the “arch” region in the substrate binding site. 45,46 According to these studies, a partial removal of the helical lid by truncation of the helix B can reduce the substrate affinity by up to 10-fold, largely by adversely affecting the off-rates. Mutagenesis experiments with different lid constructs in the DnaK chaperone have

further

highlighted the role of lid helices in controlling the stability and function of the Dnak-substrate complexes.47-49

Recent experimental studies of allosteric regulation in the DnaK

have

shown that some mutant variants, including I438A, V440A, L484A, D481A, D481K, could interfere

with distinct steps of allosteric communication and produce direction-specific

communication pathways between the NBD and SBD regions.50 Collectively, these functional studies have suggested a dynamic mechanism of substrate binding and release in which the nucleotide binding signal may be transmitted through a network of allosteric

interactions

connecting the NBD-SBD interface with the substrate binding loops.

ACS Paragon Plus Environment

5

Journal of Chemical Information and Modeling

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 6 of 114

While DnaK chaperone is by far the most extensively studied member of the Hsp70 family, structural studies of other Hsp70 proteins, such as human protein Hsp70 (HSPA1A), have only recently gained a momentum due to their therapeutic importance and growing efforts to develop effective modulators of chaperone activities.

The crystal structure of the HSPA1A complex

with a peptide substrate has shown considerable similarities to the DnaK-substrate complex, yet revealing conformational differences in the loop Lα, β that bridges the SBD-α and the SBDβ subdomains, and changes in the inter-domain linker (sometimes referred as loop LL,1) .51 Deregulations of signal transduction pathways are often linked to the Hsp70-regulated processes, and selective inhibition of the Hsp70 chaperone oncogenic pathways,

could lead to the disruption of

while achieving tumor cell specificity.52,53

multiple

These revelations have

strengthened the arguments about viability of Hsp70 as an attractive drug target, fostering the development of Hsp70 inhibitors for the treatment of cancer, particularly to combat tumor cell resistance. Despite a great interest in the therapeutic potential of Hsp70 chaperones52,53 discovery of effective Hsp70 modulators has been challenging, and only a small number of inhibitors has been identified and characterized at structural and biophysical levels.54-58 The early screening studies have discovered a series of functionalized dihydropyrimidines that could stimulate the ATPase activity of Hsp70, whereas other structurally related molecules inhibit this function.59-61 High-throughput screening (HTS) assays have identified small molecules that can adjust the folding activity of the DnaK and its co-chaperones DnaJ and GrpE.62 Structure-based design inhibitors

of adenosine-derived inhibitors has led to the first structurally characterized Hsp70 that target the ATPase binding domain.63 The

identification and molecular

characterization of VER-155008 inhibitor64 in the complex with the NBD of human Hsp70 has shown that this ligand can arrest the NBD of human Hsp70 in a specific conformation and

ACS Paragon Plus Environment

6

Page 7 of 114

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Information and Modeling

disrupt

the allosteric control.65 A small molecule, 2-phenylethynesulfonamide (PES), was

shown to bind to the C-terminal domain and inhibit the protein folding activity of Hsp70 and DnaK chaperones.65-67 Chemical biology approaches combined with structural modeling have identified a new site in the NBD of human Hsp70 chaperone68 which led to the development of several series

of allosteric

and reversible Hsp70 inhibitors.69,70

irreversible

Recently,

structural and biophysical studies have identified a new small molecule inhibitor (often referred as PET-16) that can inhibit Hsp70 and DnaK chaperones by targeting a novel C-terminal allosteric pocket.71,72

Cellular-based functional approaches have also identified and

characterized a natural product, novolactone that can affect the cellular chaperone functions through blockage of an ATP-induced substrate release and inhibition of refolding activities.73 Biochemical and structural analyses of novolactone have shown that this inhibitor is covalently bound with HSPA1A in an allosteric pocket,

which may block the inter-domain interface and

impede allosteric transitions between open and closed chaperone forms.73

High-throughput

screening of inhibitors targeting protein-protein interactions in DnaK complexes with nonenzyme partners

DnaJ and

GrpE

have identified molecules with distinct inhibitory

mechanisms.74 Computational studies have characterized conformational dynamics of the DnaK states and have offered

several atomistic models of allosteric

communications in these molecular

chaperones.75-83 Although structure and function of the Hsp70 chaperones have been actively investigated in recent years and have produced significant insights, mechanisms of allosteric regulation of chaperone activities by small molecules remained fairly elusive

due to limited

structural and dynamic information about inhibitor binding with the Hsp70 proteins.

Among

ACS Paragon Plus Environment

7

Journal of Chemical Information and Modeling

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

important questions is

to

clarify

Page 8 of 114

differences in the ligand-induced modulation of allosteric

regulation by peptide substrates and small molecule inhibitors. In the current study,

to probe allosteric inhibition of DnaK and HSPA1A chaperones, we

combined molecular dynamics (MD) simulations and binding free energy computations of the chaperone

structures with the modeling of allosteric communications and structure-based

network analysis. For network-based analysis, we employed a graph-based representation of protein structures84,85 that provides a convenient description of the residue interaction networks using topological parameters such as residue centrality.86-91 Using these approaches, we investigated and systematically compared conformational dynamics, allosteric propensities of protein residues

and organization of the residue interaction networks in

the DnaK

and

HSPA1A structures, including the apo forms and complexes with the peptide substrates and allosteric inhibitors PET-16 and novolactone (Figure 1).

To characterize

ligand-induced

changes in the allosteric interaction networks, we also determined structural organization and dynamic evolution of local stable communities in the DnaK and HSPA1A structures. show that

We

small molecule inhibitors may target central mediating sites in the residue

interaction networks and differentially modulate allosteric regulation by inhibiting specific allosteric interactions and impairing the global network connectivity. By using a mechanisticbased

analysis of distance fluctuation

maps and

residue-based modeling of allosteric

propensities, we highlight principal differences in the allosteric communication networks of the PET-16 and novolactone complexes. The results of this study unveil molecular details of allosteric inhibition in the Hsp70 proteins from both atomistic and network-centric perspectives, offering a plausible mechanism that could rationalize distinct functional effects of PET-16 and novolactone on the substrate binding and chaperone cycle.

ACS Paragon Plus Environment

8

Page 9 of 114

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Information and Modeling

Materials and Methods Structure Preparation 500 ns MD simulations were

performed independently for the DnaK-SBD and HSPA1A-SBD

crystal structures including the unliganded structures, the complexes with the peptide substrate NRLLLTG, and the complexes with the allosteric inhibitors PET-16 and novolactone respectively (Figure 1). In addition, 5 independent 200 ns simulations were carried out for each of the studied systems to ensure computations of the system properties by averaging over multiple and

independent trajectories and also to minimize potential biases that may be

encountered in a single trajectory analysis.

We simulated

all reported apo structures of

DnaK-SBD, including the 1.4 Å crystal structure of the DnaK-SBD (residues 389-607, pub id 4F01),92 the NMR structure of the DnaK-SBD-β (residues 393-507, pub id 1DG4),24 and the three independently solved crystal structures of the unliganded DnaK-SBD (2.36 Å DnaKSBD-A, pub id 4R5J, 1.75 Å DnaK-SBD-B, pub id 4R5K, and 2.97 Å DnaK-SBD-C, pub id 4R5L).71 We report the results of MD simulations for the apo-DnaK structure (residues 381-605, pub id 4R5K). MD simulations of the DnaK-SBD complexes included DnaK-NRLLLTG X-ray structures (pub id 1DKZ, 1DKY),18

the 1.98 Å DnaK-SBD complex with the NRLLLTG

peptide ( residues 384-605, pub id 4R5I),71 and the DnaK-SBD complex with the PET-16 inhibitor (residues 378-595, pub id 4R5G).71 MD simulations were also conducted for the HSPA1A crystal structures in the unliganded form (residues 395-543, pub id 4WV5),73,the complex with the peptide substrate NRLLLTG (residues 386-613, pub id 4PO2),51 and the complex with novolactone (residues 395-543, pub id 4WV7).73 In simulations we employed several

SBD constructs with different SBD-α helical lid length,

in the presence and in the

absence of the intermolecular linker. The crystal structure of the DnaK-SBD complex with the

ACS Paragon Plus Environment

9

Journal of Chemical Information and Modeling

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

substrate (residues 384-605)

Page 10 of 114

provided coordinates for the complete SBD construct, with all

helices in the helical lid structurally resolved. In contrast, in the crystal structure of the DnaKSBD complex with PET-16 (residues 378-595), helices A and B of the helical lid do not form a continuous helix and the coordinates of the helices C, D, an E were missing. Structural changes, induced by lid truncations may have important biological consequences,

as the DnaK variant

which lacks the helical bundle (helices C-E of the lid) revealed a substantially less stable secondary structure for the helix B, causing the reduced substrate-activated ATPase rate (about two times lower than that of full-length DnaK).25 All crystal structures were obtained from the Protein Data Bank (RCSB PDB www.rcsb.org).93 Structure preparation process consisted of several stages and included checking and correcting for errors, adding missing residues and atoms, filling valences with hydrogens, predicting pK values for titratable amino acids, assigning predefined partial charges and radii to all atoms, and generating force field parameter and topology files for MD simulations. Here, we described the main components and details of structure analysis and preparation. First, protein residues in all studied crystal structures were systematically inspected for missing residues and

protons.

Hydrogen atoms and missing residues were initially added and assigned according to the WHATIF program web interface.94,95 Protonation states for titratable groups were initially determined using WHAT IF pKa calculations in aqueous solution at pH 7.0.96,97 The initial assignments of protonation states were checked for consistency and subjected to a second round of refinement using H++ approach and web server (http://biophysics.cs.vt.edu/H++).98,99 Using this approach, missing hydrogens were added to the studied protein structures based on predicted pKa of titratable groups. These calculations are based on the continuum solvent approach100

in the framework of the Poisson- Boltzmann (PB) model.97 The unresolved

ACS Paragon Plus Environment

10

Page 11 of 114

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Information and Modeling

structural segments in the crystal structures of the DnaK complex with PET-16 (residues 596610) , the apo form of HSPA1A (residues 544-613) and the HSPA1A complex with novolactone (residues 544-613) were modeled using MODELLER program.101 modeling and a template-based loop

By combining homology

prediction approaches ModLoop102,103

ArchPRED104 we reconstructed and refined

and the

missing loops in the DnaK and HSPA1A

complexes. The chaperone structures were then optimized using the 3Drefine method that is based on an atomic-level energy minimization using a composite physics and knowledge-based force fields.105 The crystallographic coordinates of the DnaK-PET16 and HSPA1A-novolactone complexes were used to extract the inhibitor conformation, followed by optimization of the inhibitor geometry with density functional theory at the B3LYP/6-31G level using the Gaussian 09 package.106 In this stage of structure preparation, we followed the strategy of comprehensive and compatible parameterization for the protein and inhibitor atoms CHARMM General Force Field (CGenFF)-compatible parameters.107

utilizing a set of The atomistic charge

parameters for PET-16 and novolactone inhibitors were initially probed using ParamChem approach and webserver by assigning parameters based on analogy with a set of CGenFF parameters.108,109

The inhibitor parameters were subsequently generated an refined using the

restrained electrostatic

potential (RESP) charge fitting procedure110-112 on the electrostatic

potentials produced by single-point quantum mechanical calculations at the Hartree-Fock level with a 6-31G* basis set. The adopted protocol is consistent with the accepted best practices in the employment of CHARMM force field for simulations of protein-inhibitor complexes that include RESP fitting procedure or extracting charges directly from quantum mechanical calculations. A VMD plugin ffTK was used to complete computation of partial atomic

ACS Paragon Plus Environment

11

Journal of Chemical Information and Modeling

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 12 of 114

charges and parameters for small molecule inhibitors PET-16 and novolactone.113 This protocol is based on the CGenFF force field and parameterization strategy that allows for compatibility between

CHARMM22 force field parameters for proteins

topology and parameter

file information

for PET-16 and

Information). It is worth noting that a set of

and CGenFF parameters (see novolactone

in Supporting

conceptually similar tools such as

ANTECHAMBER114 and ACPYPE115 could simplify the automatic generation of topology and parameters and assignment of partial charges for small molecules from quantum mechanical computations that are compatible with the General Amber Force Field (GAFF).116 MD Simulations 500 ns MD simulations were carried out using NAMD 2.6 package117 with the CHARMM22 force field118,119 and the explicit TIP3P water

model.120

The employed MD protocol is

consistent with the overall setup that was described in details in our earlier studies.81,121,122 The initial structures were solvated in a TIP3 water box with the buffering distance of 10 Å from the protein surface to the solvation boundary using Solvate and Autoionize Plug-ins in VMD program.123 Sodium and chloride ions were added to neutralize the simulated systems and to ensure an ion concentration of 150 mM ionic strength. Solvated structures were minimized for 50,000 steps using fixed heavy atoms followed by relaxation for another 50,000 steps without applying any restraints. The following protocol was used before the production stage of MD simulations. All atoms of the complex were first restrained at their crystal structure positions with a force constant of 10 Kcal mol-1 Å-2. The system was subjected to the following simulation annealing to ensure the proper equilibration. The temperature was increased from 0K to 500K at a rate of 1K per 1ps and was kept at 500K for 500ps. The temperature was then decreased from 500 K to 300K at a rate of 1K per 1ps and was kept at 300K for additional 500ps.

ACS Paragon Plus Environment

12

Page 13 of 114

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Information and Modeling

500 ns MD simulations were run in an NPT ensemble with a time step of 2 fs keeping the temperature at 300 K through a Langevin thermostat with a damping coefficient of 0.5 ps-1. The constant pressure was maintained at 1 atm with a Langevin-piston barostat.124 Long-range nonbonded van der Waals interactions were treated using an atom-based cutoff of 12Å with switching van der Waals potential beginning at 10Å. The smooth particle mesh Ewald (PME) method125,126 was employed to treat the long-range electrostatics. Numerical integration was performed using the leap-frog Verlet algorithm with 2fs time step. The SHAKE algorithm was employed to constrain the covalent bond lengths involving hydrogen atoms in the solute and solvent molecules.127 Principal component analysis (PCA) of the MD conformational ensembles was used to obtain the dynamic cross-correlation matrix between residues. PCA of simulation trajectories was based on the set of backbone heavy atoms (N, Cα, Cβ, C, O) to determine the essential dynamics of the chaperone structures. The calculations were performed using the CARMA program128 and PCA_NEST web service.129 Coarse-Grained Modeling of Functional Dynamics: The Gaussian Network Model Coarse-grained approaches and elastic network models (ENM) such as Gaussian network model130-132 (GNM) combined with the normal mode analysis133,134 can efficiently probe functional movements by reducing nodes where the

protein structure

to a network of uniformly connected

native interactions in the equilibrium structure determine conformational

dynamics of the system.135-137 In the GNM representation, the protein structure is modeled as a network of N nodes identified by the α-carbon atoms of proteins where the fluctuations of each node are isotropic and Gaussian. For nodes i and j , equilibrium position vectors, are R 0i and R 0j , equilibrium distance vector is R 0ij , instantaneous fluctuation vectors are ΔR i and ΔR j , and

ACS Paragon Plus Environment

13

Journal of Chemical Information and Modeling

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 14 of 114

instantaneous distance vector is R ij . The topology of the protein structure is described by N × N Kirchhoff matrix of inter-residue contacts, Γ, where the off-diagonal elements are −1 if the nodes are within a cutoff distance and zero otherwise. Bonded and nonbonded pairs of residues located within an interaction cutoff distance rc =7.0 Å are assumed to be connected by springs with a uniform spring constant γ, which is the single parameter (force constant) of the Hookean potential.130-132 The equilibrium dynamics of the structure results from the superposition of N-1 nonzero modes found by the eigenvalue decomposition of Γ. By assigning a spring constant, γ, to all contacts, the cross-correlations between the fluctuations ΔR i and ΔR j of residues (nodes) i and j are computed. The GNM normal modes are found by diagonalization of the Kirchhoff matrix, Γ = UΛUT .

ΔR i  ΔR j 

3kBT



[Γ-1 ]ij 

3kBT



[UΛUT ]ij 

3k BT



N 1

λ k 1

-1 k

[ukukT ]ij (1)

Here, U is a unitary matrix, UT  U1 of the eigenvectors u k of Γ and Λ is the diagonal matrix of eigenvalues λ k .The elements of the kth eigenvector u k describe the displacements of the residues along the kth mode coordinate, and the kth eigenvalue, λ k , scales with the frequency of the kth mode, where 1 ≤ k ≤ N – 1. KB is the Boltzmann constant, T is the absolute temperature. The normal mode analysis generated eigenvalues λ k and eigenvectors u k that were obtained using the WEBnm@ approach.138 The normal mode analysis included computation of 10 low frequency modes using an approximation that was proven to be highly efficient and accurate in computations of low-frequency motions.139 The mean square fluctuations of a given residue can

ACS Paragon Plus Environment

14

Page 15 of 114

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Information and Modeling

be then evaluated as a sum over the contributions of all modes. The fluctuation of the ith atomic degree of freedom along the eigenvector u k reflects the mobility of residue i in the kth mode.

(ΔR i )2 

3kBT



N 1

 (λ k 1

u ukT )ii (2)

-1 k k

The residue-based conformational mobility profiles in the essential space of low frequency modes were estimated using oGNM computation of structural dynamics based on the GNM approximation.132 Using this approach, we assessed cumulative contributions of the three slowest modes to the essential residue mobility and computed the mean square displacements of protein residues driven by these slow modes. Protein Stability Calculations A systematic alanine scanning of protein residues in the DnaK and HSPA1A structures was performed using FoldX approach.140,141 We utilized a graphical user interface for the FoldX calculations142 that was implemented as a plugin for the YASARA molecular graphics suite.143 If a free energy change between a mutant and the wild type (WT) proteins ΔΔG= ΔG (MT)-ΔG (WT) > 0, the mutation is destabilizing, while when ΔΔG