Studies on the Molecular Mechanism between HDAC8 and Inhibitory in Different Bioactivities by Molecular Docking and MD Simulations①

LIANG Zhen YAN Wen-Li LI Hong-Mei LI Ying ZHANG Rong

(School of Pharmacy,Guangdong Pharmaceutical Uniνersity,Guangzhou 510006,China)

ABSTRACT HDAC8 is an important target for the treatment of many cancers and other diseases.To develop potent and selective HDAC8 inhibitors,molecular docking and molecular dynamics (MD) simulations were employed for investigation of the mechanism of HDAC8 inhibitions containing hydroxamic acid group.Compound 1 with high activity and compound 2 with low activity were selected for comparative study.Compound 1 formed a stronger chelation with Zn ion and was more stable in the HDAC8 pocket than compound 2.Residues HIS-180,ASP-178,ASP-267,and GLY-140 played a critical role in securing the position of compound 1.Both the head and tail of compound 1 formed strong hydrogen bonds with ASP-178,facilitating the ZBG of compound 1 close to the Zn ion so that they formed permanent chelation during the simulation period.The Cap group of the compounds with branch and long chains was advantageous to form interaction with active pocket opening.What's more,based on the results of this study,three innovative recommendations for the design of highly active HDAC8 inhibitors were presented,which will be useful for the development of new HDAC8 inhibitors.

Keywords:HDAC8 inhibitors,molecular docking,molecular dynamics (MD) simulations,inhibition mechanism;DOI:10.14102/j.cnki.0254-5861.2011-3126

1 INTRODUCTION

Histone deacetylases (HDACs) are ‘eraser’ enzymes that cleave acetate from acetylated e-amino groups of lysine in histones and other proteins[1-3].In humans,there are 18 kinds of HDACs divided into class I,II,III and IV[4].The HDAC8 belongs to class I and is found in both nucleus and cytoplasm.It possesses 377 amino acids (aa,42 kDa)[5-8]and is a zinc-dependent HDAC[9,10].

HDAC8 is mainly studied as a cancer target and is highly expressed in a variety of cancer tissues compared to normal tissues[11-14].Inhibiting HDAC8 is effective in the treatment of cancers containing T-cell lymphoma[15],Myeloid leukemia[16],Colon cancer[17],Breast cancer[18],Gastric carci-noma[19],Neuroblastoma[20],Hepatocellular carcinoma[21],Urothelial cancer[22],etc.Besides,inhibition of HDAC8 can fight the parasitic infections[23,24],viral infections[25,26]and Cornelia DE Lange syndrome (CdLS)[27].

Till date,only six HDAC inhibitors have successfully passed clinical trials.However,the inhibitors interact with all classes of HDACs and lead to adverse effects.Therefore,it is important to develop specific HDAC8 inhibitors.Most compounds containing hydroxamic acid group have a good inhibitory effect on HDAC8.Ingham et al.[28]have found a highly potent HDAC8 inhibitor namely OJI-1 (HDAC8 IC50=0.8 nM),a kind of hydroxamic acid.It was identified as the most selective and potent HDAC8 inhibitor till now over all other HDACs (HDAC1 to HDAC7 and HDAC9).

Currently,no specific HDAC8 inhibitor has been marketed.HDAC8 inhibitors known have either low activity or high toxicity.Therefore,it is still of great significance to continue to develop HDAC8 specific inhibitors with high activity and low toxic and side effects on human body.In addition,the specific inhibitory mechanism of inhibitors containing hydroxamic acid group on HDAC8 is still not clear,which needs further study.

Nowadays,computational aid drug design (CADD) is a successful technique in drug discovery[29,30].Recently,Zhou et al.have used molecular docking and MD simulation to study the mechanism by which three series of novel urushiol derivatives interact with HDAC8[31].Uba et al.have reached a conclusion that 1T64 (HDAC8 in complex with TSA) and similar crystal structures may be good candidates for HDAC8 structural dynamics studies and inhibitor design[32].

In this work,a high active compound and a low active compound were selected for molecular docking and MD simulations.The results of the two compounds were compared to reveal the mystery of inhibition of HDAC8 by the compounds containing hydroxamic acid group.

2 MATERIALS AND METHOD

2.1 Materials

In this work,the highest active compound 1 and the lowest active compound 2 were selected among 61 HDAC8 inhibitors containing hydroxamic acid group and synthesized by Zhang et al.[33-35](The structures are shown in Fig.1).Compound 1 (IC50:0.063uM) is 23.5 times more active than the SAHA (Vorinostat,IC50:1.48uM[36]) and 193 times more active than compound 2.The receptor is the crystal structure of HDAC8 (1T64) downloaded from Protein Data Bank (https://www.rcsb.org/).

Fig.1.Structure of SAHA (Vorinostat)[35].Compound 1 with high activity and compound 2 with low activity.The ZBG (Zinc ion Binding Group),linker and Cap group (surface recognition domain) are indicated in green,red and blue,respectively.These three parts are the HDAC8 pharmacophore[9]

2.2 Molecular docking

AutoDock Vina 1.1.2 package software was performed for molecular docking.The receptor is the crystal structure of HDAC8,and the ligand is compounds 1 and 2 in the database.The PDBQT format for the receptor and ligand molecules was prepared in MGL tools 1.5.6.The maps of binding site were calculated with 15 × 15 × 15 grid points of 1 ? spacing.35.The molecular docking was repeated 5 times to get 100 conformations of the ligands.The ligands with optimal conformation in them were selected as the initial conformation of the molecular dynamics simulation[37].

2.3 Molecular dynamic simulation

The molecular dynamics simulation was carried out using Gromacs5.1.2 software.The GROMACS 96 force field is used to calculate proteins and compounds.The system was placed into a periodic cube with conditions of 298 K and 101.3 kPa and solvated with simple point charge (SPC) water molecules.Sodium ions were added to the box to satisfy the electrical neutrality condition.The calculation is regulated and conserved by isothermal-isochoric ensemble (NVT) and isothermal-isobaric ensemble (NPT).Once the equilibrium calculation of 100 ps was completed under the NVT and NPT ensembles,a 50 ns molecular MD simulation was performed[38].

3 RESULTS AND DISCUSSION

3.1 Docking analysis

HDAC8 active site consists of a narrow and long hydrophobic tunnel,and a Zn ion is in the end of the hydrophobic tunnel[39].The docking results of compounds 1 and 2 with HDAC8 are presented in Fig.2,providing the mechanism about the compounds inhibiting HDAC8.The optimal docking score of compounds 1 and 2 was -8.8 and -7.8 kJ/mol,respectively.This indicated that the structure of compound 1 was more favorable for binding to the active site than that of compound 2.The compounds both chelated with Zn ion in active pocket,which was crucial mechanism of compounds containing hydroxamic acid group[9].The ZBG (Zn ion Binding Group),Linker and Cap groups of compound 1 (Fig.1) interacted with the residues of the active site.The Cap group of compound 1 was branched and interacted with both sides of the pocket opening,so compound 1 was embedded in HDAC8 like an iron nail,which fixed its position and made it be longer chelation with Zn ion.This is clearly and intuitively shown in Fig.3a.However,only compound 2’s ZBG and Linker interacted with the receptor.The Cap group had no interaction with the pocket opening due to its lack of branching and short length.Thereby,it is like an embroidery needle and easy to get out of the pocket (Fig.3b).

Fig.2.(a) Interactions between the active site of HDAC8 (color:green) and compound 1 (color:wheat).Metal chelating bond,H-bonds,π-π T-shaped and π-anion are depicted as cyan,yellow,red and blue dotted lines,respectively.(b) The Zn ion is shown as dark grey sphere.Interactions between the active site of HDAC8 (color:green) and compound 2 (color:wheat).Metal chelating bonds,H-bonds,π-π T-shaped and π-sulfur are depicted as cyan,yellow,red and orange dotted lines,respectively.The Zn ion is shown as dark grey sphere

Fig.3.(a) Conformation of compound 1 bound to the active site.(b) Conformation of compound 2 bound to the active site

The hydrophobic interaction of the compounds with active site is shown in Fig.4.The compounds had hydrophobic interactions with many residues.The ZBG of the compounds formed a lot of hydrogen bonds (H-bonds),while the Linker and the Cap group formed lots of hydrophobic interactions.Moreover,significantly more residues formed hydrophobic bonds with compound 1 than 2,suggesting that hydrophobic interaction was also the reason why compound 1 binds more strongly to HDAC8 than 2.

Fig.4.Hydrophobic interaction of complex system (a) HDAC8-compound 1 and (b) HDAC8-compound 2

3.2 Analysis of MD simulations

3.2.1 Binding mode of complex system

The root mean square deviation (RMSD) was applied for investigating the stability of complexes and conformational change.As shown in Fig.5,RMSD value of HDAC8 with compounds 1 and 2 were lower than the RMSD of HDAC8 without ligand,which showed that the two compounds stabilized HDAC8 and formed more stable complexes with HDAC8.RMSD of HDAC8 with compound 1 had little fluctuation and HDAC8 with compound 2 exhibits larger fluctuation after the preliminary stability (25～50 ns),which manifested that the interaction between compound 1 and HDAC8 was stronger than that between compound 2 and HDAC8.

Fig.5.RMSDs of HDAC8

The binding free energy of the two compounds to the receptor is shown in Fig.6.The Van der Waals energy (VdW energy) accounted for a greater proportion of the binding free energy between compound 1 and HDAC8.However,the electrostatic energy (Elec energy) made up a higher proportion of the binding free energy between compound 2 and HDAC8.The Vdw energy of compound 1 was three times that of the Elec energy,and the Elec energy of compound 2 was slightly greater than the VdW energy.In addition,the binding energy of compound 1 to the receptor was stable within 50 ns with little fluctuation,whereas it is opposite for 2.What's more,the binding energy of compound 1 was generally lower than that of 2.Therefore,1 has stronger interactions with HDAC8.The conclusion is consistent with that of RMSD.

Fig.6.Fluctuation of binding free energy of complex (a) HDAC8-compound 1 and (b) HDAC8-compound 2

Binding energy contribution values indicated the impor-tance of residues interacting with the ligand.Generally,the smaller it is,the greater the contribution to the combination is.It is bad for the combination if it is positive.Binding energy contribution values of the Zn ion and 9 residues are in Table 1.Binding energy contributions of Zn ion was -63.6 kJ/mol in the interaction between compound 1 and HDAC8,and it is 48.9 kJ/mol in the interaction between compound 2 and HDAC8,which revealed that the Zn ion had strong interactions with compound 1.On the contrary,it is repulsive to compound 2,which is bad for it to combine with compound 2.Chelation of the inhibitor with Zn ion was vital to inhibit the activity of HDAC8,so compound 1 was much more active than 2.

Table 1.Binding Energy Contributions between HDAC8 Residues and Compounds 1 and 2

3.2.2 Interaction between ligands and HDAC8

HIS-180,ASP-178 and ASP-267 chelated with Zn ion in HDAC8 to fix the position of Zn ion[40].Interestingly,the O atom of the hydroxamic acid group of compound 1 not only chelated with Zn ion but also formed strong H-bonds with the three residues (Fig.7).On the other side relative to the O atom,a H atom formed a constant H-bond with the O atom on GLY-140.At the cap group of compound 1,the O atom of the imide group had strong H-bonds with the H atom of imidazole ring of HIS-180.All the H-bonds did not break during the entire simulation period,indicating that the force was firm.It was the indestructible H-bonds that bonded to HDAC8 in all directions like iron chains,binding the ligand and receptor tightly together.In addition,compound 1 also formed a temporary H-bond with PHE-208,GLN-263,GLY-303,etc.,further fixing the position of compound 1 in the active pocket.Therefore,compound 1 strongly inhibits the catalytic activity of HDAC8.

Compound 2 interacted with HDAC8 rather than compound 1 with HDAC8 (Fig.8).During the simulation period,TYR-306,His180 and GLY-151 residues formed H-bonds with compound 2 for a long time and played a role in fixing the position of the hydroxamic acid group of compound 2.However,GLY-151 and TYR-306 far from the Zn ion resulted in the O atom on the carbonyl group of hydroxamic acid group also far away from it,so it was difficult to bond with the Zn ion.And the Cap group of compound 2 did not interact with receptor and was easy to move,which further damaged the stability of the ZBG.Therefore,compound 2 was so unstable that it was easy to leave the active pocket,resulting in its low activity.

3.2.3 Interaction with Zn ion

Fig.7.Variation of HDAC8-compound 1 in the active sites during the MD trajectory

Fig.8.Variation of HDAC8-compound 2 in the active sites during the MD trajectory

The literatures showed that the main mechanism of HDAC8 inhibitor containing hydroxamic acid group is the chelation with Zn ion[15,41,42].During the simulation period,only can the O atom on the carbonyl group of hydroxamic acid group formed chelation with the Zn ion.Therefore,we figured out that the distance between the two atoms varies within 50 ns (Fig.9).Firstly,the distance between the O atom of compound 1 and Zn ion maintained a steady level and small values (around 0.2 nm) in the 50 ns.This indicated that the interaction between compound 1 and the Zn ion was too firm to damage,while the distance between the O atom of compound 2 and Zn ion was greater than 0.4 nm.Such distance was difficult to chelate with Zn ion.What's more,this distance was larger in the later simulation period.

Fig.9.Distance between Zn and the O atom on the carbonyl group of hydroxamic acid group

3.2.4 Analysis of H-bonds

H-bonds played an important role in the interaction between the ligands and receptor.In general,the more H-bonds,the more stable the complex.The probability of the compounds forming H-bonds with the HDAC8 is shown in Fig.10.There were at most six H-bonds formed between compound 1 and the receptor,while compound 2 had only five.Compound 1 had the highest frequency of two and three H-bonds with the receptor,while one and two H-bonds for 2,indicating compound 1 was more likely to form more H-bonds with HDAC8 than compound 2,that is,the complex of HDAC8-compound 1 is more stable than HDAC8-com-pound 2.

Fig.10.Relative frequency of the number of hydrogen bonds formed in complex (a) HDAC8-compound 1 and (b) HDAC8-compound 2 within 50 ns

3.2.5 Analysis of the covariance matrix

The covariance matrix was applied as a supplement for understanding the variation of residues in HDAC8.The matrix value of the white region was 0.The value of blue region was negative,representing negatively correlated motion.The value of red region was positive,showing the positively correlated variation.The covariance matrix of HDAC8 C-alpha is displayed in Fig.11.The matrix of compound 1 has more red regions than 2.In addition,the maximum value of red region of HDAC8 with compound 1 was 0.274 and 0.211 with compound 2.These indicated that 1 prompted more positively correlated variation in HDAC8.

Fig.11.Covariance matrix of (a) HDAC8 with compound 1 and (b) HDAC8 with compound 2

3.3 Molecule design of new HDAC8 inhibitors

In general,the compounds containing hydroxamic acid group are mostly inhibitors of HDAC8,but their activity varies greatly due to different structures.According to the research results,a few conclusions may improve the activity of this class of compounds.First,the Cap group of the compounds with branch and long chains is advantageous to form more interactions with opening the active pocket,by which the compounds embedded in the HDAC8 like an iron nail fix their positions (as shown in Fig.3).Second,when the ZBG and Zn ion are in the chelating positions,an O,N or F atom should be added to the part of compounds near the HIS-180 in order to form strong H-bond with the H atom of imidazole ring of HIS-180.For example,the O atom of the imide group of compound 1 forms H-bond with the H atom.This design ensures ZBG to remain in the end of the tunnel and keep chelation with Zn ion for a long time.Third,the branch chains of the Cap group can be designed to benzene rings,so that they formπ-π,π-anion,π-sulfur and other interactions with the residues in pocket opening to further stabilize the position of the compounds.

4 CONCLUSION

The inhibition of HDAC8 activity by compounds containing hydroxamic acid group is mainly dependent on chelation with Zn ion.In order to maintain the interaction for a long time,the compounds stabilized in the active pocket to make sure the O atom of the hydroxamic acid group was close to the Zn ion.For the receptor,HIS-180,ASP-178,ASP-267 and GLY-140 played a significant role in stabilizing the ligand position.For ligand,the O atom of the imide group of the Cap group forming H-bond with the H atom of imidazole ring of HIS-180 played a crucial role in fixing the ligand's position.The Cap group of compound 1 is a branching structure composed of two long chains,which facilitates its interaction with the opening of active pocket and stabilizes the position of the tail of compound 1.It is like a nail in the active pocket.However,compound 2 has a short tail and no branching structure,which is unfavorable for stabilizing its location in the active pocket.Furthermore,several conclusions about the design of new high-activity HDAC8 inhibitors have been drawn:the cap group of the compounds should be designed as a long branching structure containing benzene rings.And an O,N or F atom should be added to the part of compounds near HIS-180 when the ZBG and Zn ion are in the chelating position.These conclusions have significance implications for the development of highly active HDAC8 inhibitors in the future.