Molecular Dynamics Simulations Reveal the Mechanisms of Allosteric Activation of Hsp90 by Designed Ligands

Controlling biochemical pathways through chemically designed modulators may provide novel opportunities to develop therapeutic drugs and chemical tools. The underlying challenge is to design new molecular entities able to act as allosteric chemical switches that selectively turn on/off functions by modulating the conformational dynamics of their target protein. We examine the origins of the stimulation of ATPase and closure kinetics in the molecular chaperone Hsp90 by allosteric modulators through atomistic molecular dynamics (MD) simulations and analysis of protein-ligand interactions. In particular, we focus on the cross-talk between allosteric ligands and protein conformations and its effect on the dynamic properties of the chaperone’s active state. We examine the impact of different allosteric modulators on the stability, structural and internal dynamics properties of Hsp90 closed state. A critical aspect of this study is the development of a quantitative model that correlates Hsp90 activation to the presence of a certain compound, making use of information on the dynamic adaptation of protein conformations to the presence of the ligand, which allows to capture conformational states relevant in the activation process. We discuss the implications of considering the conformational dialogue between allosteric ligands and protein conformations for the design of new functional modulators.

Protein functions are determined by their internal dynamics and are fine-tuned by the interactions with effectors of different chemical origins [1][2][3] . Most chemical interventions on proteins and enzymes focus on the identification of small molecule inhibitors designed to abrogate their activities 4 , or on the development of agonists that induce active-like conformations of the target and subsequent signaling responses by binding to known orthosteric sites 5 .
In contrast, the rational search for allosteric activators of protein functions is still relatively unexplored 6,7 . However, allosteric activators may provide new opportunities to develop chemical tools to probe the role of protein activity in cellular phenotypes: indeed new information can be gained by turning on a signaling pathway starting from a specific node of the underlying protein network 7,8 . Allosteric activators can in principle modulate receptor function while still allowing the possibility of orthosteric agonist or antagonist binding. In this context, allosteric ligands can act as affinity modulators by changing the affinity of the orthosteric ligand for the receptor, as well as induce conformational changes that in turn modulate the efficacy of the orthosteric ligand in determining cellular responses.
Moreover, activators can help identify and characterize allosteric sites and mechanisms for the discovery of novel drug candidates: indeed, since allostery represents one of the most relevant means to regulate protein function, drug discovery can be extended to target regulatory protein pockets. It has previously been shown that sequence variation in allosteric pockets could aid the design of highly specific drugs that would bind only specific members of a protein family. Finally, allosteric ligands can potentially modulate the balance between the target protein conformations that are presented to other interacting proteins, generating the possibility to perturb protein-protein interactions (PPIs) through low molecular weight ligands 9,10 .
Defining the principles of ligand recognition and allostery is thus important for understanding the molecular links between small-molecule binding and observed effects. This would in turn facilitate the generation of novel chemical tools and therapeutic agents. Experimental approaches based on X-ray crystallography, tethering, enzymatic, biochemical and phenotypical assays have been used to characterize the activation of enzymatic functions determined by small molecules 6,7,9,10 . Despite these exciting advances, experimental techniques are still limited in their potential to provide insight at an atomic level into the allosteric activation process itself. To understand activation mechanisms at atomistic detail, we can turn to theoretical and computational approaches 10,11 . The latter methods can help address open questions such as the role of variations in protein structures and internal dynamics incurred by binding at an allosteric site, and the possibility to quantitatively correlate ligand binding, dynamic-changes and observed modulation of a certain activity 12 .
In this paper, we will address the problem of rationalizing the observed allosteric activation of the 90 kDa Heat Shock Protein (Hsp90) molecular chaperone through a set of newly designed ligands 12 . Members of the Hsp90 family are hub proteins that control the crossroads of fundamental pathways required for cell development and maintenance. The disregulation of Hsp90 functions, often associated with its overexpression, has been shown to be at the basis of disease states, such as cancer and neurodegeneration [13][14][15][16][17][18][19] . At the molecular level, the activities of Hsp90 are regulated by ATP binding and hydrolysis: the type of bound nucleotide selects the conformational sub-states and the transitions among them, which ultimately determine the chaperone's functional properties 20 .
From a structural point of view, Hsp90 is a homodimer (see schematic representation in Fig. 1), whereby each protomer consists of three distinct structural domains: an N-terminal regulatory Domain (NTD), where the ATP binding site is located, a Middle Domain (M-domain) composed of a Large (LMD) and a Small (SMD) subdomain, which completes the ATPase site necessary for ATP hydrolysis and binds client proteins, and a C-terminal Domain (CTD) that is involved in dimerization [21][22][23][24][25] . A bipartite client binding site made of residues from both the MD and CTD has been identified 26,27 . This structural organization is common to the various Hsp90 homologues from different organisms and cellular compartments [21][22][23][24][25]28 .
Although the exact mechanisms of coupling between ATP-binding/hydrolysis and the chaperone cycle remain partly elusive, combining available crystal structure information and solution state data, Agard and collaborators have developed a general model that recapitulates the main traits of nucleotide dependent conformational changes in Hsp90. In particular, ATP binding in both NTDs from the two protomers stabilizes a closed NTD-dimerized activated state, characterized by an excess structural strain that is ultimately responsible for the population of an asymmetric conformation 28,29 . Hydrolysis of one of the two ATPs could then cause a change in symmetry (either to a symmetric state or possibly an opposite asymmetric state), providing a mechanism to directly couple hydrolysis to rearrangement of the client binding residues identified in the MD:CTD region [21][22][23]28 . A simple schematics of the mechanism is reported in Fig. 1A.
With the aim of shedding light on the roles of the ATPase in regulating the conformational and cell-functional properties of the chaperone, we recently designed a library of small molecules that activate Hsp90 ATPase by targeting an allosteric site located 65 Å away from the active site 12 . The latter, located at the MD:CTD border in a region largely overlapping with the client binding site at the interface between the two protomers, was identified through the characterization of long-range dynamic coordination patterns based on the analysis of residue-pair distance fluctuations [30][31][32] . The characterization of the stereochemical and dynamic properties of the allosteric site constituted the basis for the optimization of the allosteric ligands, generating compounds able to accelerate the chaperone enzymatic activity up to 6 times.
Herein, we examine the origin of the allosteric activation of Hsp90 and set out to lay the basis for the development of a predictive model of ligand-induced activation of the protein's functional dynamics. Through the use of comparative analyses of molecular dynamics (MD) simulations of different complexes of the chaperone in the presence of allosteric activators, and known inhibitors as controls, we provide molecular-level insights into how ligand binding at an allosteric site can affect protein structure, dynamics and, consequently, enzymatic activity. We also provide a first, dynamics-based, quantitative model that links the structure of the allosteric ligands to experimentally measured activation effects. Finally, we discuss the implications of our work for activator design.

Results and Discussion
To start building a structural and dynamic model of the mechanisms of Hsp90 activation induced by the allosteric ligands, we carried out MD simulations of closed ATP-bound Hsp90 in complex with compounds 1, 4, 10, 12, 16, 18 ( Fig. 1B) at the allosteric site (Fig. 1A, top-right), and compared the results with those obtained in the absence of compounds (ATP-only). The structures, ATPase and closure rate stimulation effects for the compounds are reported in Table 1. We focused on the closed conformation of Hsp90 as a model of the activated state of the protein. In this context, our aim is not full sampling of conformational changes or the exhaustive exploration of the highly complex free energy landscape of the protein and its modulation by the ligands. Our goal is to shed light on the microscopic details of Hsp90 internal dynamics induced by specific ligands that can ultimately be linked to the activation of a specific functional state, namely the one related to ATPase activation. This is expected to determine dynamic signatures for defined ligand states of the chaperone, providing information on allosteric mechanisms.

Structural analysis: Allosteric Activators Stabilize Closed Compact States of the Dimer.
We first analyzed the frequency distribution of distances between the centers of mass of the N-terminal domains of the two protomers in the presence of a certain allosteric ligand, compared to the ATP-only bound states. Only backbone atoms were considered in calculation. This parameter informs on the propensity of Hsp90 to populate the closed, N-terminal dimerized state in the presence of the allosteric ligand. For the complexes of Hsp90 with ATP and allosteric activators, the peaks of the frequency distributions are at similar or slightly shorter distances than the ATP-only simulation, indicating an enhanced propensity to populate closed conformations ( Fig. 2A). This latter is considered to represent the activated state of the chaperone 31,32 .
In contrast, control simulations run with known inhibitors (Novobiocin, DHN1, DHN2) show that the peaks for the inhibitors are consistently centered at larger NTD-NTD distances (Fig. 2B), indicating a propensity for stabilizing a more open, inactive state 31,32 .
Overall, within the limitations of MD sampling, these results indicate that the allosteric activators tend to maintain the protein in an ensemble of structures similar to the ones observed for the ATP-only case, where the N-terminal organization remains tightly closed. The population of the latter, initiated by ATP (or ATP analogs), is a hallmark of Hsp90 activation and represents the catalytically active state of the protein that allows ATPase activity. In contrast, inhibitors appear to favor less compact conformations by shifting the peak distributions by about 0.5 nanometers (5 Å) towards more open structures, which may be associated to the selection of non-reactive states 31,32 . As a caveat, it must be stated that these considerations clearly have a mostly qualitative value and more efficient sampling methods, based on enhanced sampling techniques, coupled to longer simulation timescales might be needed to quantitatively profile the effects of allosteric ligands on Hsp90 free energy landscape.

Allosteric Activators Distort Hsp90
Inducing an Asymmetric State. We next analyzed the large-scale distortions determined by the allosteric ligands. To do this, we evaluated the time evolution of the angle formed by the principal axis of the N-M domains of each protomer in simulated snapshots with the  principal axis of the same N-M domain in the 2CG9.pdb closed structure, derived from yeast and used a reference, after alignment of the structures on the CTDs. In the absence of the allosteric ligands, the two protomers overall show a similar behavior with a limited distortion relative to the starting conformation. In contrast, the presence of the activators consistently determines a clear distinction between the two protomers, with protomer B, which contacts the allosteric ligand more extensively, being significantly more distorted than protomer A (Fig. 3). The only exception to this behavior is represented by the complex formed with compound 16, whose profile is highly similar to that of the ATP-only complex. Interestingly, 16 is the derivative showing the least effect on ATPase activity (see Table 1). The rationale for this different behavior can be found in the specific structure of the carbohydrate moiety (β -L-Fucose) in 16: the cis configuration of the two hydroxy groups at positions C-3 and C-4 and the lack of an OH group at position C-6 of the sugar do not favor the formation of stable hydrogen bonding interactions with the side chains of E477 and D503 on one of two protomers, which are located on opposite sides of the sugar ring. These residues can be ideally targeted by trans OH groups on atoms C-3 and C-4, as well as by on OH substitution on C-6, as in the case of compounds 1, 4, 10, and 12 (see Fig. 4). In the case of compound 18, the positively charged tertiary amino-group acts as a bridge between the two negatively charged side chains of E477 and D503, through electrostatic interactions. Such network of interactions induces strain in the 3D structure of the protein, which eventually results in the observed distortion in one of the two protomers, namely the one to which the two interacting E477 and D503 side chains belong. The structural details at the basis of the activity of the different ligands are also described in our previous paper on allosteric design 12 .
In this framework, the distortions induced by the activators could favor the population of the asymmetric state proposed by Agard and coworkers as a high-energy intermediate for the activation of ATP hydrolysis 28,29 . Effects of Allosteric Ligands on Hsp90 Internal dynamics. Next, we turned our attention to the characterization of the effects of allosteric ligands on the internal dynamics of Hsp90, compared to the ATP-only case. We first characterized the global rigidity/flexibility patterns in the complexes by calculating the fluctuations of pairwise amino acid distances (this analysis was also named coordination propensity (CP) analysis in previous papers [30][31][32]. To obtain a compact representation of the data for the different compounds, we calculated the difference matrix between all residue pair fluctuations determined for the ATP-only case and all residue pair fluctuations determined in the presence of the allosteric ligands. Positive values in the matrix (black color) indicate higher flexibility in the ATP-only state, while negative ones (blue color) indicate higher flexibility in the allosteric complex. This calculation is then intended to pinpoint the distribution of rigidity/flexibility blocks within the protein. As a consequence of ligand binding, amino acids belonging to different domains may change their long-distance coordination, coupled to a modulation of the inter-domain rigid dynamics.
In the presence of the activators (Fig. 5 Upper Panel), the patterns of flexibility/rigidity are highly similar to the ones of the ATP-only case, as shown by the large white regions, the main differences being that in the allosteric complexes the NTDs from one protomer (represented as a light blue domain in the figure) are decoupled from the block formed by the Small portion of the Middle Domain (Small-Middle Domain, SMD) and the CTD (represented as a light red and green domains in the figure) of the other, as indicated by the larger blue-colored blocks. Interestingly, the split of the Middle domain into a Large and a Small subdomain is consistent with the domain separation observed in the X-ray structure of the mitochondrial Hsp90 homolog TRAP1 and with evolutionary studies indicating that Hsp90 and other GHKL family members significantly diverge at this subunit 33 . These results once more suggest a peculiar role for the SMD-CTD unit in regulating functionally oriented aspects of Hsp90 conformational dynamics.
By contrast, the high degree of internal coordination between the NTD and MD within each protomer was maintained upon compound addition, indicated by the presence of extended white blocks. As observed previously, intraprotomer coordination of the two domains that make up the active site favors structural preorganization for catalysis. In parallel, interprotomer flexibility can be aptly exploited to speed up the search for the closed active state, consistent with the observed increases in ATPase.
In contrast, control simulations run in the presence of coumarin-based inhibitors indicate a uniform increase in the overall protein flexibility, likely favoring distortions of both protomers and detachment of the two NTDs ( Fig. 5 Lower Panel). This enhanced flexibility upon addition of a compound is quite unusual as compound binding normally stabilizes proteins, reducing dynamics. Thus, it may be a hallmark of allosteric inhibitors to work by increasing dynamics, favoring transitions away from the ground and catalytically active states. In the case of The design and selection of small molecule inhibitors targeting well-defined enzymatic catalytic sites is now routinely accessed through virtual screening. On the other hand, experimental screening approaches can be used to select for a certain protein activity. The Wells lab demonstrated the possibility to select for specific kinase and caspase activators by combining a tethering approach with crystal structure studies 7,34 . Work by the Buchner lab proved the feasibility of using high-throughput FRET based assays for the discovery of artificial activators of Hsp90 35 .
However, the computational design of allosteric protein activators still represents a challenging task, in general and for Hsp90 in particular. The dynamic nature of the allosteric binding sites, combined with the fact that in many cases no direct correlations between the ligand binding affinity values and any effects on chaperone functions or cellular activity can be observed, makes the problem of structure-based activators design particularly complex [36][37][38] .
To progress along this difficult but fascinating avenue, we set out to explicitly consider the dynamic cross-talk between the chaperone and allosteric activators as an integral part of the development of a structure-dynamics-activity relationship (SDAR) model. In other words, we made direct use of the dynamic adaptation of protein conformations to the presence of the ligand (and vice versa) in deriving a simple model relating the physico-chemical properties of the small molecules to their activating effects. On this basis, we chose to describe binding in terms of an ensemble of ligand poses in complex with an ensemble of protein structures: in this way, protein and ligand conformations that could play a role in the activation process may be captured and recapitulated in a simple model.
To achieve this goal, we adopted a Docking over Multiple Receptor Structures (DMRS) approach, originating from the ATP-only dynamics: the 10 representatives of the most populated conformational clusters obtained from the ATP-only simulation were extracted from the MD trajectory and compounds 1, 4, 10, 12, 16, 18 were re-docked into each of the 10 representative structures (see Materials and Methods) (Fig. 6 exemplifies this for  (Table 1 -Avg. LE on ATP-only column) strongly correlates for each activator with the experimental ATPase rate (Table 1 -ATPase), thus representing a predictive descriptor for the evaluation of new modulators. The combination of DMRS and LE computation provides a convenient common structural framework to design, dock and evaluate novel compounds, without the need for additional computationally expensive MD simulations. As a control, for each activator, we also evaluated the average LE from dynamic trajectories of the respective Hsp90 complexes, using the representatives of the most populated structural clusters for each dynamics (Table 1 -Avg. LE on Clusters) and on random snapshots extracted from the equilibrated parts of the trajectories (Table 1 -Avg. LE on Snapshots). The resulting LE values for each ligand were then averaged and correlated to the respective experimentally determined ATPase stimulations (Table 1). Interestingly, the correlation coefficients resulted to be higher (− 0.90) using the docked poses in the 10 protein structures from the ATP-only simulation, than using the representative structures of the most populated conformational clusters of the complex dynamics (− 0.70), or using the representative snapshots extracted from the trajectories of the same dynamics (− 0.57). These data indicate that LE evaluated on multiple protein structures can aptly capture the salient determinants of observed ATPase activation profiles and allow to extract quantitative structure-dynamics-activity relationships for a series of modulators. Control calculations using only one single best pose yielded poor correlations in the order of − 0.10. This result underlines the importance of using multiple poses of the allosteric ligands in contact with different conformations of the chaperone to facilitate the generation of a model that quantitatively correlates the structures and dynamics of ligand/protein complexes to experimentally observed activation of ATPase activity (Fig. 6). To the best of our knowledge, this is one of the first fully rational structure-dynamics-activity models for the modulation of chaperone functions induced by allosteric small molecule ligands.
We also evaluated the correlations between average LE values and FRET closure rates. A significant correlation of − 0.76 could be obtained only when using the poses from DMRS (values of − 0.40 and − 0.30 were obtained for the structures obtained from the complex dynamics clusters and snapshots, respectively). We speculate that in this latter case, larger conformational changes are relevant to correctly describe the closure process and sampling limitations play a key role.
Though based on a limited number of compounds, which is a necessary condition given the high computational costs of all-atom MD simulations, our data overall support a model whereby the allosteric small molecules stabilize specific protein conformations selected from the ensemble around the closed activated state, possibly providing a path towards precise control of Hsp90 reactivity.
Indeed, the analysis of the effects of compounds on the structural distortions of Hsp90 highlighted relevant patterns for activation. Consistent with the characterization of the time evolution of the angle formed by the principal axes of the N-M domains of each protomer compared to the 2CG9.pdb crystal structure, an evident asymmetry in the two protomers emerges in the presence of the allosteric activators, with protomer B (the one in contact with the variable functional groups installed on the benzofuran scaffold) consistently showing a larger distortion than protomer A (Fig. 7).
In this context, it is possible to hypothesize that different ligand poses modulate the ensemble of protein conformations responsible for enzymatic activities. The coupling of dynamic ligand binding to distinct active protein conformations defines thus the mechanism for allosteric activation.

The Cross-Talk between Protein and Allosteric Ligand Determines Allosteric Activities.
The possibility of characterizing the dynamic properties of a protein's active state can provide a rational basis for designing chemical stimulators of its function. The main challenges in this framework involve the advancement of drug design to explicitly include the characterization of protein functional dynamics, the development and experimental testing of compounds able to directly modulate the functional substates of the target through selected functional groups, and the investigation of the mechanistic consequences of allosteric binding on the functionally oriented dynamic properties of the protein. Meeting these challenges can directly contribute to advancing the role of chemical design in the investigation of complex biomolecular systems. The long-term perspective is in fact the definition of methods for the rational development of new chemical entities able to modulate and not only inhibit protein activities linked to cellular functions.
Our work demonstrates that the cross-talk between allosteric ligands and protein conformations can dynamically modulate protein activities. In the conformational selection model of binding and activity, a small molecule can preferentially target and stabilize key functional sub-states from an ensemble of conformations available on the protein free energy landscape. This concept was for instance at the basis of the design of Gleevec, a potent inhibitor of kinase activities, which binds differently to the active conformation of a kinase (Syk) than to the inactive conformation of another (Abl) 40 .
In this context, we have previously combined conformational selection concepts and studies of nucleotide-regulated dynamics for the discovery of an allosteric site in Hsp90 [30][31][32] . Next, we used the knowledge of the stereoelectronic properties of the allosteric site to select a small molecule lead binder. This was then evolved into a small library of derivatives showing the ability to allosterically activate Hsp90 ATPase and closure dynamics. More potent activators were then designed through the rational modification of the functionalities decorating the benzofuran scaffold to exploit potential productive interactions with the network of charged amino acids in the allosteric pocket, resulting in compounds able to accelerate Hsp90 enzymatic activity and closure kinetics by up to 6 times 12 .
Building on these results, herein we set out to generate a comprehensive structural, dynamic and activity model to link the presence of specific allosteric ligands to the experimentally observed activation of Hsp90 functional stimulation. Binding of allosteric activators determines the stabilization of the closed N-terminal dimerized In parallel, analysis of Hsp90 internal dynamics using the distance fluctuation Coordination Propensity (CP) parameter highlights lower interprotomer coordination in the presence of the allosteric activators. In this case, an asymmetric behavior emerges as well: in protomer A, in fact, high coordination entails all the NTD and M-large regions. The subdomains of protomer B, on the other hand, show patterns similar to those observed for the ATP-only cases.
Although never previously observed with a eukaryotic cytosolic Hsp90, the proposed formation of an asymmetric, activated state is consistent with recent results by the Agard Lab based on the crystal structures of the mitochondrial Hsp90 Trap1 and studies on the bacterial Hsp90 HtpG. These results indicate that after ATP binding at the NTD, a high-energy asymmetric state is populated which favors sequential ATP hydrolysis by the two protomers 28,29 . The ligand-induced increased flexibility coupled to the differential structural and dynamic profiles between the two protomers observed here can favor the global rearrangements that lead to and stabilize the formation of the high-energy asymmetric catalytically competent state, eventually resulting in ATPase stimulation.

Conclusions
Overall, our findings establish that the interactions between the allosteric ligands and the closed active conformations of Hsp90 can direct functional dynamics towards the preferential population of compact asymmetric states, in which conformers of the chaperone that are optimally poised for ATPase are selected, collectively determining the observed biological activity. This model is consistent with the description of allosteric activation events recently proposed by Nussinov and coworkers 41 , whereby different states of a protein can be separated by a certain surmountable free energy barrier, and binding of the allosteric effectors shifts the population eliciting allosteric effects. In this context, the small molecule effector acts as a conformational catalyst that favors the protein's transition to the active state, lowering the barrier for the population shift. Moreover, our work also highlights the functional importance of transient asymmetric states in otherwise symmetric protein assemblies.
While still hampered by the inherent limitations of MD sampling for such a large and complex system as Hsp90, our work begins to define how allosteric ligands can elicit modulated responses in their target proteins, while at the same time shedding light on the underlying mechanisms of protein dynamics and setting the stage for the use of computational approaches to rationally design chaperone activators. These aspects can have important implications in the development of allosteric modulators with graded activities for a number of relevant targets, including for instance nuclear receptors, kinases, caspases and GPCRs. In all these contexts, the availability of chemical tools that allow the control and modulation of the functional dynamics of proteins (and their networks) can have fundamental impact on physiology studies and the development of novel mechanism-based therapeutics.

Materials and Methods
Molecular Dynamics Simulations. The starting structures of MD simulations of Hsp90 in complex with the allosteric ligand compounds were obtained from docking calculations, using the most representative protein conformation of previous MD simulations of Hsp90 in the ATP bound state 30 . The receptor preparation was performed using the "Protein Preparation Wizard" suite of the MAESTRO program (see www.Schrodinger.com). Bond orders and atomic charges were assigned, and hydrogen atoms were added. The assignments of protonation states for basic and acidic residues were based on the optimization of hydrogen bonding patterns. The final minimization of the protein was performed with the default parameters of the Preparation Wizard. The allosteric binding site previously defined 30  The simulation and the analysis of the trajectory were performed using the GROMACS 45 software package using the GROMOS96 53A6 force field 46 . In order to remove any bad contacts, the complex was initially minimized in vacuo by multiple minimizations (200 steps steepest descent plus 200 steps conjugate gradient). Then the system was solvated in a tetrahedral box large enough to contain the protein and 1.0 nm of solvent around it on each side. The SPC water model was used for solvation 47 , and Na+ counterions were added to ensure electroneutrality. In order to allow the solvent molecules to relax around the solute, the system was minimized with position restraints on the protein and ligand and just minimizing the positions of water and ions (500 steps steepest descent plus 200 steps conjugate gradient). After this minimization step, the whole system was equilibrated by a 100 ps of MD run in the NVT ensemble to allow the temperature reaching a plateau at 300 K 48 . The first equilibration run was followed by a second 100 ps run in the NPT ensemble to adjust the pressure and the density of the system, by weak coupling to a bath of constant pressure (P0 = 1 bar, coupling time τ P = 0.5 ps) 48 . The first 20 ns of the trajectory were not used in the subsequent analysis in order to minimize convergence artifacts. Equilibration of the trajectory was checked by monitoring the equilibration of quantities, such as the root-mean-square deviation (rmsd) with respect to the initial structure, the internal protein energy, and fluctuations were calculated on different time intervals. The electrostatic term was described by using the particle mesh Ewald algorithm 49 , the LINCS algorithm 50 was used to constrain all bond lengths and the SETTLE algorithm 51 was used for water molecules. A dielectric permittivity, ε = 1, and a time step of 2 fs were used. All atoms were given an initial velocity obtained from a Maxwellian distribution at the desired initial temperature of 300 K. In all simulations, the temperature was maintained close to the intended values by weak coupling to an external temperature bath with a coupling constant of 0.1 ps. The proteins and the rest of the system were coupled separately to the temperature bath 48 . Each Hsp90/allosteric molecule complex was run for 100 ns. Three replicas were run for each complex.
All analyses were carried out with the tools available in the Gromacs suite and with in house developed software. For the contact analysis reported in Fig. 2, a contact is considered to be present if the distance between any pairs of atoms from two residues belonging to different protomers was lower than 6 Å.
Protein internal dynamics. To analyze the impact of each compound on protein internal dynamics, we made use of the previously introduced distance fluctuation analysis. For each MD trajectory in complex with an allosteric ligand (activator or inhibitor) we computed on the time interval 25-100 ns (the first 25 ns trajectory are removed in order to avoid equilibration artifacts), the matrix of distance fluctuations, corresponding to the CP parameter 31 : where d ij is the (time-dependent) distance of the Cα atoms of amino acids i and j and the brackets indicate the time-average over the trajectory. Notice that CP is invariant under translations and rotations of the molecules and, unlike the covariance matrix, does not depend on the choice of a particular protein reference structure. The CP matrix, and various quantities derived from it, can be used to characterize the salient elasticity and plasticity properties of a protein undergoing structural fluctuations. The presence of specific coordination patterns and quasi-rigid domains in the protein of interest should reflect in specific properties of CP. In fact, pairs of amino acids belonging to the same quasi-rigid domain should be associated with much smaller distance fluctuations than amino acid pairs in different domains.
The CP values obtained from three different 100 ns MD trajectories for each allosteric ligand were averaged. The same was done with the matrix corresponding to the ATP-only case. We first characterized the global effects of ligands on chaperone dynamics by calculating the fluctuations of pairwise amino acid distances on all possible residue pairs: the data obtained for the simulations of Hsp90 in complex with both ATP plus allosteric ligands were subtracted with the ones obtained for the chaperone in complex with ATP only, and reported in the matrices analyzed in the text. This allws to identify residues that change their coordination upon ligand binding.
Axes and angles. In order to investigate the internal rearrangement of the protein induced by the compounds, in particular the asymmetry of the protein and the distortion within each monomer, we defined two different axes of inertia of every protomer and then we determined the two angles between them. To do this, all the trajectories were sampled every 500 ps starting at 50 ns, obtaining 100 structures per simulation. Using the program Chimera (v 1.9) 52 each monomer of all protein conformations was superimposed on the CTD of the corresponding monomer in the crystal structure of Hsp90 (PDB code: 2CG9). Two axes per monomer were defined: the first along the entire protomer and the second comprising only the NTD and Middle domains. To determine the axes, in all cases we used only residues forming the conserved secondary elements along the simulation time.
We then measured all the angles formed by NM axis of this protomer and all the other protomers, considering the average on the three replicas. Finally we calculated the difference (NM monomer B -NM monomer A) of each structure, reporting the asymmetrical rearrangement of the protein along the simulation time due to the presence of different compounds. In order to highlight the internal deformation of each protomer, we also measured the angle between the entire monomer and its NM axis considering, as stated above, the average values of the replicas and then the differences (chain B -chain A).

Structural Cluster Analysis.
Structural cluster analysis of the trajectories was carried out by the method described by Daura and coworkers 53 . Briefly: count the number of neighbors using a cut-off of 0.2 nm for RMSD between the optimal backbone superimposition of the sampled structures, then take the structure with the largest number of neighbors with all its neighbors as cluster and eliminate it from the pool of clusters. This procedure is repeated for the remaining structures in the pool. The RMSD was calculated on the backbone atoms of the C-terminal domain. We kept the first ten structures for the docking calculations, representing more than 95% of the structure variability.
Quantitative Structure-Dynamics-Activity Relationships (QSDAR). The first ten most populated structures resulting from the cluster analysis of the MD simulations of Hsp90 in the ATP bound state described in 30 without any ligand at CTD were used for docking calculations in the DMRS approach. To extract control data from the complex dynamics, the representative structures of the 6 most populated clusters resulting from the combined analysis of the MD trajectories produced for the closed ATP-bound Hsp90 in complex with compounds 1, 4, 10, 12, 16, 18 (Table 1 -clusters) or the structures of 6 snapshots (for each ligand complex) extracted from the equilibrated part of the same MD trajectories (Table 1 -snapshots) were used. In DMRS, re-docking calculations were carried out using the grid described in the Molecular Dynamics Simulations paragraph. Docking calculations were carried out with the program Glide 43 in Extra Precision mode (XP) with standard OPLS-AA (2001) force field and, in order to refine the conformational search, 20000 poses per ligand for the initial phase of docking were saved with a scoring window for keeping initial poses within 200 kcal/mol, keeping up to the best 800 poses per ligand for successive energy minimization.
In the second and third case, LE was simply calculated as a Glide score given calculated with MAESTRO on the respective structures.
At the end of calculation the experimental data of the ATP-ase activity and FRET measurements were correlated by Pearson coefficient with the Average Ligand Efficiency calculated based on the Glide-score, given by Glide in all ten conformations of the complexes (results are shown in Table 1).