The heat shock protein 90 (HSP90) is a chaperone protein responsible for maturation and activity of a variety of key client proteins involved in cellular processes   . This protein is overexpressed in tumour cells, and causes uncontrolled proliferation of transformed cells   . HSP90 function is dependent on its ability to bind and hydrolyse ATP at the N-terminal domain   . An intrinsic ATPase activity is required for the operation of a functional chaperone cycle leading to the stabilization of client proteins   .
The investigation of HSP90 as anticancer target has arisen from evident progress in clinical trials using the natural inhibitors geldanamicyn (GMD) and radicicol as well as the potent inhibitor 17-allylaminogeldanamycin (17-GMD) . However, several toxicities, extensive metabolism and cell resistance hampered the clinical progression. There is, therefore, a growing interest in design of new therapeutic inhibitors to overcoming these limitations    . This spurred development of new classes of pharmacophores using varied pharmacological approaches. Computer-assisted drug development is currently the most used approach as it requires less time and capital. Computer-assisted drug development based on molecular docking, molecular dynamic simulation and quantitative structure activity relationships (QSARs) have been used to predict the antagonist activity of various heterocyclic compounds against the target protein HSP90  . Investigations on compounds containing selenium are of some interest as it has been previously reported that selenium reduced the risk of numerous cancers like prostate, colon, live, lung and mammary cancers .
The present work explored the binding ability of GMD selenoderivatives at the N-terminal domain of human HSP90 using molecular docking. HSP90 was derived from the Protein Data Bank (PDB code: 1YET). The docking was performed with the GOLD 5.6 program (Genetic Optimization for Ligand Docking) implemented in the Cambridge Structural Database (CSD). Gold fitness scores and interaction energies were calculated to determine the affinity and the thermodynamic stability of the binding, respectively.
2. Materials and Methods
2.1. Preparation of HSP90
The crystal structures of human HSP90 N-terminal domain with the co-crystallized inhibitor GMD were retrieved from the Protein Data Bank (PDB codes: 1YET, 1.9 Å) and imported into Accelrys Discovery Studio 2019 and visualized with Hermes 1.6 implemented in the Cambridge Structural Database (CSD). Survey of the inhibitor binding site was carried out with the reference of amino acid residues of the ATP-N-terminal domain as previously reported . Hydrogen atoms were added to the protein for correct ionization and tautomeric states of amino acid residues. Co-crystallized water molecules, forming one or more hydrogen bonds to either the protein or the ligand in the protein binding site, were kept in the calculation, whereas the remaining water molecules and inhibitors were extracted before the docking. HSP90 was prepared for docking without minimization of protein energy. The molecular docking was performed after conversion of the HSP90 using GOLD 5.6 program implemented in CSD.
2.2. Preparation of Ligands
Geldanamycin (GMD) was the co-crystallized ligand in the inhibitor site of HSP90. The selenoderivatives of GMD, labelled as GMDSe (i) were designed using Mercury 3.10 program implemented in CSD by replacing an oxygen atom (O) according to Scheme 1 by a selenium atom (Se) in geldanamycin at the position given by the i label .
Scheme 1. Geldanamycin (GMD) et its GMDSe (i) derivatives (i = 1, 2, ∙∙∙, n).
The ligand conformations were generated and optimized after structure conversion from 2D to 3D as illustrated by the GMDSe1, GMDSe6, GMDSe7and GMDSe8 examples in Figure 1 using Mercury 3.10. Only energetically least conformers were selected for further use.
Figure 1. The 3D structures of the GMDSe1, GMDSe6, GMDSe7 and GMDSe8 selenoderivatives.
The molecular electrostatic potentials (MEP) of ligands were calculated at the Hartree-Fock theory (HF) implemented in the Gaussian 09 program  using the basis sets 6-31G.
2.3. Molecular Ligand-Protein Docking
The GOLD 5.6 genetic algorithm was used to generate the bioactive binding poses of GMD selenoderivatives in the active site of protein HSP90. Before docking HSP90 with ligands, the reproduction of experimental poses using the N-terminal domain of 1YET and 49 additional HSP90 complexes with co-crystallized ligands derived from PDB (Table 1) was performed to validate the program GOLD 5.6. A good reproduction was observed, if the value of Root Mean Square Deviation (RMSD) was lower or equal to 2 Å .
Table 1. PDB codes and resolutions (R) of the 50 HSP90 complexes used for validation of GOLD 5.6.
The docking of HSP90 with the GMD selenoderivatives was performed using the default GOLD fitness function ChemPLP by considering the flexibility of the ligand and specific amino acid residues of the protein site. Greater the ChemPLP fitness score, better the binding affinity. Rescoring with Chemscore was performed in order to determine the free enthalpy of ligand binding (ΔGbind) according to the equation below (Scheme 2).
Scheme 2. The equation describes the free enthalpy of binding . The ΔG0 (−5.48 kcal∙mol−1) and S (bar) components are experimental and entropic terms of the binding energy, while ΔGHB, ΔGMet, ΔGLip, ΔGRot represent energy terms related to hydrogen-bonding, metal, lipophilic and rotation interactions, respectively.
Each simulation was performed ten times yielding ten docked conformations unless three of the ten poses were within 1.5 Å RMSD of each other. Conformer pose with the lowest energy was considered as the binding conformation in the protein site. The program Accelrys Discovery Studio 2019 was used to model non-bonded polar and hydrophobic contacts in the inhibitor site of HSP90.
3. Results and Discussion
3.1. Validation GOLD 5.6 Potential for Ligand-Protein Docking
The distribution of RMSD values derived from the docking of 50 HSP90 examples with co-crystallized ligands was performed to validate use of the GOLD 5.6 program. In addition, RMSD values calculated with GOLD 5.6 for the docking of 1YET, 2BMS and 2VCI were compared to those previously reported by Lauria et al.  using the program IFD (Induced Fit Docking). The superposition of experimental poses of the co-crystallized ligands and those derived from the molecular docking with GOLD 5.6 was also carried out in order to validate the potential of GOLD 5.6 to modelling ligand-protein interactions.
The distribution of RMSD values determined for 50 examples of HSP90 complexes is shown in Figure 2.
Figure 2. RMSD distribution derived from docking of 50 HSP90 examples with co-crystallized inhibitors.
As it can be seen from the Figure 2, experimental poses of co-crystallized inhibitors were reproduced by the docking with GOLD 5.6 for 44 HSP90 examples (88%) with RMSD values less than 2 Å, which is recommended by the literature as a limit value for a good pose reproduction      .
Table 2 reports the RMSD values from docking of 1YET, 2BMS and 2VCI using GOLD 5.6 program and those reported by Lauria et al.  using the program IFD (Induced Fit Docking).
Table 2. RMSD values for 1YET, 2BMS and 2VCI complexes.
The results of docking of 1YET, 2BMS and 2VCI complexes listed in Table 2 reveal that the values of RMSD calculated with GOLD 5.6 were smaller than those with IFD  suggesting GOLD 5.6 as a powerful program for ligand-protein docking.
Figure 3 illustrates the experimental pose of co-crystallized GMD from PDB (1YET) and typical superposition of the pose derived from the molecular docking with GOLD 5.6.
Figure 3. Complex of HSP90 with GMD: (a) HSP 90 complex with co-crystallized inhibitor derived from PDB (1YET); (b) Superposition of the experimental pose of GMD (in blue) and the pose from molecular docking (in green) with GOLD 5.6:RMSD = 0.40 Å.
The experimental poses and the poses derived from the docking were superposed with RMD values less than the 2 Å cut-off level as illustrated in Figure 3(b). This also confirms the potential of the GOLD 5.6 1 program to modelling ligand-protein interactions.
3.2. HSP90 Docking with GMD Selenoderivatives
Figure 4 shows the contact areas involved in the interaction between the N-terminal domain site of human HSP 90 and the selenoderivatives of GMD as derived from the docking with GOLD 5.6.
Figure 4. HSP90 complexes with contact areas involved in the binding of GMDSe6, GMDSe7 and GMDSe8 to the site of N terminal domain.
Figure 4 demonstrates that hydrogen-bonding and lipophilic interactions were involved in the bonding of the selenoderivatives of GMD to the inhibitor site of HSP 90. The affinity of the ligands to the binding site values was determined based on calculation of GOLD fitness scores and interaction energies. The values of GOLD fitness scores, interaction energies and RMSD and are collected in Table 3.
Table 3. GOLD fitness score, RMSD, the ΔGBind free enthalpy of binding and its ΔGHB, ΔGLip, ΔGRot and S(bar) component terms in kcal∙mol−1 (ΔG0 = −5.48 kcal∙mol−1 and ΔGMet = 0).
The results of docking in Table 3 show that the interactions were generated with values of RMSD less than 2 Å. Values of GOLD fitness scores and free enthalpies of ligand binding to the HSP90 site (ΔGBind) suggest that GMD as well as selenoderivatives interacted with high affinity to the N-terminal domain of HSP90 leading to thermodynamically stable complexes. According to the GOLD fitness scores, the binding affinity of HSP90 to the ligands decrease in the sequence GMDSe6 > GMDSe7 > GMD > GMDSe8 > GMDSe5 > GMDSe2 > GMDSe4 > GMDSe9 > GMDSe3 > GMDSe1, while the decrease sequence is GMDSe6 > GMDSe7 > GMD > GMDSe8 > GMDSe5 > GMDSe2 > GMDSe4 > GMDSe3 > GMDSe9 > GMDSe1 according to free enthalpies of binding (ΔGBind). GOLD fitness scores as well as ΔGBind suggest GMDSe6 and GMDSe7 as the best ranked ligands, while GMDSe1 as the last ranked. The values of the ΔGHB and ΔGLip component terms related to hydrogen-bonds and lipophilic interactions, respectively, reveal that these interactions were the main forces stabilizing ligand-protein complexes.
The ΔGLip component term related to lipophilic interactions seems to be enhanced for the binding of GMDSe7, GMDSe6 and GMDSe2 to HSP90 suggesting that the replacement of oxygen atoms by divalent selenium atoms increases their hydrophobic contact areas. One could also expect that covalently double-bonded selenium atoms in GMDSe7, GMDSe6, GMDSe3 and GMDSe2 interact with the O and N nucleophile atoms of amino acid residues as divalent selenium atom is a sigma-hole donor . However, sigma-hole interactions, well known as short, weakly attractive and linear contacts between regions of lower electronic density (σ-hole) of donor atoms and electron lone pairs of nucleophile atoms acting as sigma-hole acceptors , were not found in the inhibitor site of HSP90.
It is noteworthy that ΔGBind markedly decreases when replacing the oxygen atom in the O1=C or O9=C group of GMD by monovalent selenium atom in the Se1=C or Se9=C group to obtain GMDSe1 and GMDSe9, respectively. Indeed, the value of ΔGBind for the binding of GMD decreases in comparison with the values of ΔGBind for the binding of GMDSe1 and GMDSe9 from −39.65 kcal∙mol−1 to 31.91 kcal∙mol−1 (19.5%) and −35.80 kcal∙mol−1 (9.7%), respectively. Likewise, the value of the ΔGHB component term decreases from −30.19 kcal∙mol−1 for the interaction of HSP90 with GMD to −19.79 kcal∙mol−1 (34.4%) and −26.86 kcal∙mol−1 (11%) for the interactions with GMDSe1 and GMDSe9, respectively. Bibelayi et al.  also reported that monovalent selenium atom is less hydrogen-bond donor than oxygen atom. The authors demonstrated based on crystallographic data derived the CSD and quantum mechanical calculations that the ability of the monovalent Se to forming hydrogen bond stems from resonance-induced Cd+=Sed− dipoles as illustrated for Se in Scheme 3 for [C,N]C=Se, contrarily to hydrogen-bonding ability of oxygen atom, which is ascribed to Cd+=Od− permanent dipoles.
Scheme 3. Electron delocalization in (N, C) C=Se system.
The hydrogen-bonding can be ascribed to electrostatic interactions between hydrogen-bond donors and acceptors of polar amino acid (AA) residues present in the HSP90 binding site and those of GMD and its selenoderivatives, while lipophilic interactions are probably due to the contacts between hydrophobic amino acid residues of HSP9 and non-polar areas of the ligand molecules. Contact areas of hydrogen-bonding and lipophilic interactions in GMD and its selenoderivatives were evidenced by molecular electrostatic potentials (MEP) as illustrated in Figure 5.
Figure 5. Computed electrostatic potential on 0.001 a.u molecular surfaces of GMD, GMDSe6, GMDSe1 and GMDSe7.
The electrostatic potentials calculated on 0.001 molecular surfaces of GMD and its derivatives in Figure 5 feature marked positive (in blue), negative (in red) and non-polar (in grey) regions corresponding to hydrogen-bond donor, hydrogen-bond acceptor and hydrophobic areas, respectively. Figure 6 reports typical hydrogen-bonding and lipophilic contacts involved in the interaction of HSP90 with GMD.
Figure 6. HSP90 complex with co-crystallized inhibitor (GMD): (a) Complex derived from PDB (1YET); (b) Contact patterns of hydrogen-bonding and lipophilic interactions involved in the binding of the inhibitor at the N-terminal domain of HSP90.
The activity of several potential HSP90 inhibitors has been predicted based on molecular docking in conjunction with molecular dynamic simulation and quantitative structure activity relationship(QSAR) studies  . Indeed, therapeutic properties of HSP90 inhibitors have been ascribed to lipophilic and hydrogen-bonding contacts involving key amino acid residues  . The modelling of contact areas like illustrated in Figure 6 shows lipophilic and hydrogen-bonding interactions involving typical hydrophobic and polar amino acid residues like Asp93, which is crucial for ATP activity at the N-terminal domain of HSP90  . The modelling of binding interactions also showed that GMD was bonded to the hydrophobic amino acid residues like Leu48, Met98, Leu107, Phe138, Val150 and Val186 as well as the polar amino acid residues Asp93, Lys58, Lys112 and Phe138 in the inhibitor site. Interestingly, the same amino acid residues were also involved in the lipophilic and hydrogen-bonding interactions of HSP90 with all GMD selenoderivatives, except GMDSe5.
Table 4 reports hydrogen-bond donors and acceptors involved in the interactions of HSP90 with GMD and the GMDSe1, GMDSe2, GMDSe3, GMDSe4, GMDSe6, GMDSe7, GMDSe8 and GMDSe9 selenoderivatives. The interaction distances (δ) and angles (θ) were performed according to Scheme 4. The value vdW(H) = 1.10 Å determined by Rowland et al. (1996)  was used for van der Waals radius of hydrogen atom, while the values vdW(O) = 1.52Å for oxygen and vdW(Se) = 1.90 Å for selenium were taken from Bondi et al. .
Scheme 4. Geometric parameters of the hydrogen-bonding: interaction distance (𝛿) and hydrogen-bond angle (θ).
Table 4. Hydrogen-bonds derived from docking of HSP90 with GMD and selenoderivatives as ligands.
Table 4 shows that the GMDSe1, GMDSe2, GMDSe3, GMDSe4, GMDSe6, GMDSe7, GMDSe8 and GMDSe9 derivatives formed hydrogen-bonds with the key amino acid residues of HSP90 like GMD did . The values of interaction distances and angles agree well with the geometric parameters of a typical hydrogen-bond. Indeed, θ values greater than the 120˚ recommended limit . The δ bond distances are less than the rH + rA sums of the van der Waals radii of the atoms involved in the hydrogen-bonds (rH + rO =2.62 Å and rH + rSe = 3.0 Å).
Earlier clinical trials demonstrated the opportunities for using inhibition of HSP90 as target protein in the treatment of cancer. The present work was an additional support for recent advances in use of the molecular docking as a powerful tool for design and discovery of inhibitors with improved pharmacological properties. Indeed, the study revealed the ability of GMD selenoderivatives to interact with HSP90. The thermodynamic stability of the binding was mainly ascribed to hydrogen-bonds and lipophilic interactions with amino acid residues in HSP90 inhibitor site. The GMDSe6 and GMDSe7 derivatives were the best ranked GMD derivatives. They were bonded to HSP90 stronger than GMD, while maintaining lipophilic interactions and hydrogen bonds with amino acid residues like Asp93, which are catalytically crucial for therapeutic properties of HSP90 inhibitors. Inhibitors containing Se atoms may provide improved therapeutic properties as the anticancer activity of selenium has been previously reported in the literature. Therefore, a pharmacophore study of the GMDSe6 and GMDSe7 compounds should be useful to explore their therapeutic properties. Molecular dynamic simulation and QSAR studies based on ligand chemical and quantum descriptors as well as experimental investigations are suggested to assess the inhibitor activities.
Kilembe Thambwe Jason, Bibelayi Dikima Didi and Lundemba Singa Albert are grateful to the Cambridge Crystallographic Data Centre (CCDC) for funding Master and PhD fellowships, and additional financial and material support.
 Abbasi, M., Henrry, S.-A. and Amanlou, M. (2017) Prediction of New HSP90 Inhibitors Based on 3,4-Isoxazolediamide Scaffold Using QSAR Study, Molecular Docking and Molecular Dynamic Simulation. DARU Journal of Pharmaceutical Sciences, 25, 17.
 Sawai, A., Chandarlapaty, S., Greulich, H., Gonen, M., Yee, Q., Arteaga, C.L. and Solit, D.B. (2008) Inhibition of HSP90 Down-Regulates Mutant Epidermal Growth Factor Receptor (EGFR) Expression and Sensitizes EGFR Mutant Tumors to Paclitaxel. Cancer Research, 68, 589-597.
 Goetz, M.P., Toft, D.O., Ames, M.M. and Erlichman, C. (2003) The HSP90 Chaperone Complex as a Novel Target for Cancer Therapy. Annals of Oncology, 90, 1169-1176.
 Hyun, S., Woo, J.K., Yazici, Y.D., Niamh, M.O., Andrew, J.S., Trevor, T.P., Clive, W., Zisterer, D.M., Lloyd, D.G. and Meegan, M.J. (2011) Lead Identification of β-Lactam and Related Imine Inhibitors of the Molecular Chaperone Heat Shock Protein 90. Bioorganic & Medicinal Chemistry, 19, 6055-6068.
 Baby, S.T., Sharma, S., Enaganti, S. and Cherian, P.R. (2016) Molecular Docking and Pharmacophore Studies of Heterocyclic Compounds as Heat Shock Protein 90 (HSP90) Inhibitors. Bioinformation, 12, 149-155.
 Ahmad, M.S., Yasser, M.M., Sholkamy, E.N., Ali, A.M. and Mehanni, M.M. (2015) Anticancer Activity of Biostabilized Selenium Nanorods Synthesized by Streptomyces Bikiniensis Strain Ess_amA-1. International Journal of Nanomedicine, 10, 3389-3401.
 Nissink, J.W.M., Murray, C., Hartshorn, M., Verdonk, M.L., Cole, J.C. and Taylor, R. (2002) New Test Set for Validating Predictions of Protein-Ligand Interaction. Proteins, 49, 457-471.
 Eldridge, M.D., Murray, C.W., Auton, T.R., Paolini, G.V. and Mee, R.P. (1997) Empirical Scoring Functions?: I. The Development of a Fast Empirical Scoring Function to Estimate the Binding Affinity of Ligands in Receptor Complexes. Journal of Computer-Aided Molecular Design, 11, 425-445.
 Lauria, A., Ippolito, M. and Almerico, A.M. (2009) Inside the HSP90 Inhibitors Bindingmode through Induced Fit Docking. Journal of Molecular Graphics and Modelling, 27, 712-722.
 Gabb, J. and Jackson, R.M. (1997) Modelling Protein Docking Using Shape Complementarity, Electrostatics and Biochemical Information. Journal of Molecular Biology, 272, 106-120.
 Zaheer, U.H., Sobia, A.H, Reaz, U. and Jeffry, D. (2010) Benchmarking Docking and Scoring Protocol for the Identification of Potential Acetyl Cholinesterase Inhibitors. Journal of Molecular Graphics and Modelling, 28, 870-882.
 Michael, L. and Christopher, E. (2010). Comparison of Current Docking Tools for the Simulation of Inhibitor Binding by the Transmembrane Domain of the Sarco/Endoplasmic Reticulum Calcium ATPase. Biophysical Chemistry, 150, 88-97.
 Hioual, K.S., Chikhi, A., Bensegueni, A. and Merzoug, A. (2012) Comparative Data on Docking Algorithms: Keeping the Update in the Field Knowledge. The International Journal of Accounting Information Systems, 2, 2249-0868.
 Politzer, P., Murray, J.S. and Clark, T. (2013) Halogen Bonding and Other Sigma-Hole Interactions: A Perspective. Physical Chemistry Chemical Physics, 15, 1117-11189.
 Metrangolo, P., Neukirch, H., Pilati, T. and Resnati, G. (2005) Halogen Bonding Based Recognition Processes: A World Parallel to Hydrogen Bonding. Accounts of Chemical Research, 38, 386-395.
 Bibelayi, D., Lundemba, A.S., Allen, F.H., Galek, P.T.A., Pradon, J., Reilly, A.M. and Yav, Z.G. (2016) Hydrogen Bonding at C = Se Acceptors in Selenoureas, Selenoamides and Selones. Acta Crystallographica, B72, 317-325.
 Abasi, M., Sadeghi-Aliabadi, H. and Amanlou, M. (2018) 3D-QSAR, Molecular Docking, and Molecular Dynamic Simulations for Prediction of New HSP90 Inhibitors Based on Isoxazole Scaffold. Journal of Biomolecular Structure and Dynamics, 36, 1463-1478.
 Sepehri, B. and Ghavani, R. (2019). Towards the In-Silico Design of New HSP90 Inhibitors: Molecular Docking and 3D-QSAR CoMFA Studies of Tetrahydropyrido [4, 3-d] Pyrimidine Derivatives as HSP90 Inhibitors. Journal of Medicinal Chemistry, 14, 439-450.
 Rowland, R.S. and Taylor, R. (1996) Intermolecular Nonbonded Contact Distances in Organic Crystal Structures: Comparison with Distances Expected from van der Waals Radii. The Journal of Physical Chemistry, 100, 7384-7391.