Combinatorial In Silico Strategy towards Identifying Potential Hotspots during Inhibition of Structurally Identical HDAC1 and HDAC2 Enzymes for Effective Chemotherapy against Neurological Disorders

Histone deacetylases (HDACs) regulate epigenetic gene expression programs by modulating chromatin architecture and are required for neuronal development. Dysregulation of HDACs and aberrant chromatin acetylation homeostasis have been implicated in various diseases ranging from cancer to neurodegenerative disorders. Histone deacetylase inhibitors (HDACi), the small molecules interfering HDACs have shown enhanced acetylation of the genome and are gaining great attention as potent drugs for treating cancer and neurodegeneration. HDAC2 overexpression has implications in decreasing dendrite spine density, synaptic plasticity and in triggering neurodegenerative signaling. Pharmacological intervention against HDAC2 though promising also targets neuroprotective HDAC1 due to high sequence identity (94%) with former in catalytic domain, culminating in debilitating off-target effects and creating hindrance in the defined intervention. This emphasizes the need of designing HDAC2-selective inhibitors to overcome these vicious effects and for escalating the therapeutic efficacy. Here we report a top-down combinatorial in silico approach for identifying the structural variants that are substantial for interactions against HDAC1 and HDAC2 enzymes. We used extra-precision (XP)-molecular docking, Molecular Mechanics Generalized Born Surface Area (MMGBSA) for predicting affinity of inhibitors against the HDAC1 and HDAC2 enzymes. Importantly, we employed a novel in silico strategy of coupling the state-of-the-art molecular dynamics simulation (MDS) to energetically-optimized structure based pharmacophores (e-Pharmacophores) method via MDS trajectory clustering for hypothesizing the e-Pharmacophore models. Further, we performed e-Pharmacophores based virtual screening against phase database containing millions of compounds. We validated the data by performing the molecular docking and MM-GBSA studies for the selected hits among the retrieved ones. Our studies attributed inhibitor potency to the ability of forming multiple interactions and infirm potency to least interactions. Moreover, our studies delineated that a single HDAC inhibitor portrays differential features against HDAC1 and HDAC2 enzymes. The high affinity and selective HDAC2 inhibitors retrieved through e-Pharmacophores based virtual screening will play a critical role in ameliorating neurodegenerative signaling without hampering the neuroprotective isoform (HDAC1).

Introduction

Chromatin architecture plays a decisive role in transcriptional regulation which in turn is potentially modulated by the antagonistic activity of HATs and Histone deacetylases (HDACs) (Eberharter and Becker, 2002; Qi et al., 2011). While HATs favor the transcriptional activation via chromatin decondensation, HDACs promote chromatin condensation and subsequent gene silencing (Yang and Seto, 2007). The opposing activities of HATs and HDACs regulate acetylation homeostasis that plays a crucial role in governing various gene expression programs (Ropero and Esteller, 2007; Ganai, 2016a). HDACs are conjugated enzymes modulating both histone and non-histone substrates and act as corepressors in transcriptional events (Yang and Seto, 2008; Mottamal et al., 2015; Ganai, 2016b). The 18 HDACs identified in human beings till date, have been divided into four classes based on structural resemblance to yeast HDACs (Mottamal et al., 2015). Class I HDACs mainly lack shuttling ability, are ubiquitous in distribution and includes HDAC1, HDAC2, HDAC3 and HDAC8. Unlike Class I HDACs, Class II HDACs are tissue specific in distribution and possess shuttling ability (Fischle et al., 2001). This Class is further subdivided into Class IIa and Class IIb HDACs. While Class IIa covers HDAC4, 5, 7, and 9, Class IIb encompasses HDAC6 and HDAC10 (Ganai, 2017a). Class III includes mechanistically distinct HDACs termed as Sirtuins (SIRT1-SIRT7) (Albani et al., 2010; Morris et al., 2010). Class IV includes HDAC11 as the only member and in comparison to other HDACs is least studied (Seto and Yoshida, 2014). Class I, II and Class IV HDACs are Zinc (Zn²⁺) dependent and are also known as classical HDACs. Sirtuins are NAD⁺ dependent and are associated with cellular senescence and aging (Ganai, 2016c; Watroba and Szukiewicz, 2016). Among Class I HDACs, HDAC1 and HDAC2 show the highest structural identity with each other.

Histone deacetylase inhibitors (HDACi), the small molecules interfering HDACs are emerging as potent chemotherapeutic agents. Based on structural distinction these inhibitors may be hydroxamates like suberoylanilide hydroxamic acid (SAHA), Trichostatin A (TSA); benzamide derivatives like pyridin-3-ylmethyl N-[[4-[(2-aminophenyl) carbamoyl] phenyl]methyl]carbamate (MS-275 or entinostat), 4-acetamido-N-(2-aminophenyl)benzamide (CI-994); cyclic peptides like romidepsin and (3S,6R,9S,12R)-6,9-dimethyl-3-[6-(oxiran-2-yl)-6-oxohexyl]-1,4,7,10-tetrazabicyclo[10.3.0]pentadecane-2,5,8,11-tetrone (HC-toxin); short chain fatty acids encompassing sodium butyrate, phenylbutyrate and valproic acid. These inhibitors have shown promising activity against cancer and neurodegeneration. Four HDACi are approved by US-FDA for treating distinct malignancies. SAHA has been approved for treating cutaneous T-cell lymphoma (CTCL). Romidepsin followed SAHA in gaining approval and is currently used against CTCL and peripheral T-cell lymphoma (PTCL) (Ververis et al., 2013). Belinostat, the third HDAC inhibitor has been approved for relapsed/refractory PTCL (Ververis et al., 2013; Ganai, 2016a). The fourth approved inhibitor panobinostat is currently active against multiple myeloma (Ganai, 2016c). Most of the HDACi including vorinostat are pan-inhibitors targetting HDACs of different classes (Dasmahapatra et al., 2010), few like entinostat are class selective targetting isoforms of a given class (Duque-Afonso et al., 2011) and very few like tubacin are isoform selective (Lee et al., 2015; Ganai, 2016d) targetting a single HDAC.

Histone deacetylase inhibitors (HDACi) including LAQ824 (dacinostat), pyroxamide, HC-toxin are composed of three distinct components; Zinc binding group (ZBG) which chelates Zinc (Zn) ion situated deep in the active site; Linker region connecting ZBG with cap region and interacting with active site tunnel residues; Cap region which closes the active site gate and interacts with active site rim residues. This three component concept has proved spectacular in developing potent inhibitors against HDACs. For designing isoform-selective inhibitors specific modifications in these components have been exploited (Witt et al., 2009; Ganai, 2016d, 2017a).

Implications of HDAC2 in Neurodegenerative Events

The cellular imbalance between HAT and HDAC activity alters acetylation homeostasis causing transcriptional dysregulation which in turn provides impetus to neurodegenerative signaling (Saha and Pahan, 2006). Aberrant expression of HDACs has been implicated in neurologic pathologies (Ganai et al., 2016). For instance, HDAC6 overexpression has been reported to inhibit nerve growth by deacetylating tubulin protein (Rivieccio et al., 2009). Knockout studies in mouse models have revealed HDAC2 as a key regulator of associative and spatial memory (Guan et al., 2009). While HDAC2 overexpression has been reported to impair memory performance, HDAC2 knockout mice showed improved memory performance (Guan et al., 2009; Volmar and Wahlestedt, 2015). Robust improvement in associative learning has been demonstrated on selective knockout of HDAC2 (Morris et al., 2013). Upregulated HDAC2 levels have been implicated in restraining the expression of neuroplasticity genes during neurodegeneration (Gräff et al., 2012). Mounting evidences suggest that HDAC2 plays a central role in mediating cognitive impairment (Welberg, 2012). Further, the impact of elevated levels of HDAC2 on basic excitatory neurotransmission has been well demonstrated in mature neurons indicating its role in synaptic plasticity (Akhtar et al., 2009). HIV-Tat protein induced upregulation of HDAC2 culminates in down-regulation of genes meant for synaptic plasticity and neuronal function thereby triggering HIV-associated neurocognitive disorders (HAND) (Saiyed et al., 2011). These findings establish the critical role of HDAC2 in modulating synaptic plasticity and enduring changes of neural circuits culminating in negative regulation of learning and memory. Taken together, HDAC2 is a propitious and eye catching epigenetic target for tackling neurodegenerative maladies (d'Ydewalle et al., 2012; Ganai, 2017b).

Current Impediment in Targetting HDAC2 for Neurological Disorders

Pharmacological intervention with HDACi reversed the reduced synapse number and learning impairment of HDAC2-overexpressing mice (Guan et al., 2009). Glutamate excitotoxicity has been implicated in many neurodegenerative disorders and studies have shown that inhibitors active against HDAC2 and HDAC3 offer neuroprotection against such toxicity in rat models (Durham, 2012). Unlike HDAC2, HDAC1 interacts with HDRP and facilitates neuronal survival (Bardai et al., 2012). Moreover, the catalytic domains of HDAC1 and HDAC2 share high percentage identity (94%) (Figure 1) and thus therapeutic intervention against the latter also targets HDAC1 (Corpet, 1988; Ganai, 2014). This often causes debilitating off-target effects, emphasizing the desperate need of designing on-target (isoform-selective) inhibitors for HDAC2. Recently two inhibitors BRD4884 and BRD6688 demonstrating kinetic selectivity for HDAC2 vs. HDAC1 have been synthesized (Wagner et al., 2015). Another group has used a novel de novo reaction-mechanism- based inhibitor design approach toward the discovery of selective inhibitor β-hydroxymethyl chalcone against HDAC2 (Zhou et al., 2015). Taking these facts into consideration the current study used a combinatorial in silico approach including extra-precision molecular docking, molecular mechanics generalized born surface area, molecular dynamics simulation (MDS), trajectory clustering and energetically optimized structure based pharmacophore mapping for highlighting the hotspots of inhibitors in the HDAC1 and HDAC2 binding pocket. Five inhibitors belonging to three different structural groups of HDAC inhibitors were docked against HDAC1 and HDAC2 active site. These docked complexes were subjected to MMGBSA for predicting the binding affinities of docked inhibitors. The docked complexes of top scoring inhibitors LAQ824 and HC-toxin were subject to the cutting edge MDS for 5 ns. The MDS output file of docked complexes was used as input for Desmond trajectory clustering. Seven clusters were generated for each protein-ligand complex and the cluster with maximum number of frames (more stability) was considered for creating hypothesis to highlight the critical features of inhibitor inside the active site of HDAC1 and HDAC2 enzymes.

FIGURE 1

Figure 1. HDAC1 and HDAC2 share high sequence identity (94%) at the active site. The active site residues were taken from UniProt and alignment was performed by using MultAlin and cross checked by using Clustal Omega. Percent identity was calculated by Clustal Omega.

Materials and Methods

Protein Preparation and Grid Generation

Accurate starting structures are prerequisite for successful structure based modeling. The crystal structures of HDAC1 and HDAC2 (PDB ID: 4BKX and 4LY1 respectively) retrieved from Protein Data Bank (http://www.rcsb.org) (Lauffer et al., 2013; Millard et al., 2013) were prepared using the Protein Preparation Wizard of Schrödinger package (Maestro v11.0) to ensure structural correctness (Sastry et al., 2013; Ganai et al., 2015a,b). In the first step the missing hydrogen atoms were added to crystal structures and proper bond orders were assigned. Moreover, missing side chains and missing loops were filled using the Prime. All the water molecules beyond 5 Å were deleted. In the next step, the redundant protein chains and heteroatoms were deleted. As HDACs require Zinc for their catalytic function so this heteroatom was kept intact (Ganai et al., 2015b; Sinha et al., 2016; Steinbrecher et al., 2017). Moreover, the native ligand in crystal structure of HDAC2 was kept as such and was used for grid generation in the later stage. The third stage involves the refining of protein structures to make them suitable for subsequent steps. During this process, the structures are optimized and the water molecules with <3 hydrogen bonds to non-waters are deleted. This was followed by minimization in which heavy atoms were converged to Root mean square deviation (RMSD) of 0.30 Å. Grid generation was performed using the cocrystallized ligand as centroid in case of HDAC2 (Glide v7.3). However, in case of HDAC1 lacking the cocrystallized ligand, grid generation was done by specifying the residues interacting with active site Zinc (Sastry et al., 2013; Ganai et al., 2015a).

Ligand Preparation

Computational methods like molecular docking require correct 3D molecular models as initiating point. However, many compounds available in compound databases exist as 2D molecular structures and thus accurate 2D−3D conversion is a critical progenitor to computational analysis. LigPrep incorporated in Schrodinger package generates a single, factual, energy minimized 3D molecular structure with correct chiralities from the given input. Besides, it eradicates mistakes in ligands in order to enhance the accuracy of downstream events including molecular docking (Ganai et al., 2015a; Van Den Driessche and Fourches, 2017). LigPrep can optionally produce many structures from single input structure with various ionization states, stereochemistries, tautomers, and ring conformations. The structure files of 5 ligands used in the study were retrieved from PubChem with PubChem CID's provided in Supplementary Figure 1. These ligands were prepared for molecular docking using the predefined LigPrep (Maestro v11.0) (Kalyaanamoorthy and Chen, 2013; Ganai et al., 2015b). Ligands were desalted; metal binding states were generated as their receptors are Zinc dependent enzymes. Protonation states of ligands and the associated energy penalties were predicted using the Epik (Epik v3.8), a computer program based on Hammett and Taft methodology (Shivakumar et al., 2010). The preparation parameters were kept identical for all ligands to avoid any bias arising due to differential ligand preparation.

Pose Validation by Self-Docking

Among the various methods used for validating the docking programs, the pose validation method is globally used. Crystal structure of HDAC2 with its cocrystallized ligand was prepared using the protein preparation wizard followed by separation of ligand. The ligand was then redocked against the host target using extra-precision (XP) flexible docking protocol and the RMSD was calculated between crystal and redocked pose (Sandor et al., 2010).

Molecular Docking

The molecular docking was performed using the Glide (Grid-based Ligand Docking with Energetics) of Schrödinger package. The ability of Glide to identify hits for lead optimization and guide understanding of critical interactions apart from desolvation effects influencing receptor-ligand binding, has contributed markedly to its epidemic acceptance. The prepared ligands were docked against grid generated receptors HDAC1 and HDAC2 in an extra precision (XP) flexible mode (Glide v7.3) (Friesner et al., 2006). Identical docking conditions were set for both HDACs to avoid any discrepancy due to differential parameters. The GlideScore (GSore) representing affinity of ligands against receptors was obtained from pose viewer file of docked complexes (Kalyaanamoorthy and Chen, 2013; Ganai et al., 2015b).

Molecular Mechanics Generalized Born Surface Area (MMGBSA)

The Prime/MMGBSA tool is widely used for estimating the relative binding affinity of various ligands. The binding energy of ligands calculated by MMGBSA is expected to align well with the experimental binding affinity, especially of a congeneric series. Binding free energy calculations of docked complexes was performed using the Prime module of Schrodinger package. The pose viewer files generated in Glide XP docking were used as input for MMGBSA. The binding free energy was calculated in the frozen state and default dielectric constants, solvation model (VSGB) and OPLS3 force field was used (Lyne et al., 2006; Ganai et al., 2015b). Prime MMGBSA performs five fundamental energy calculations; optimized free receptor (Receptor), optimized free ligand (Ligand), Optimized complex (complex) in addition to receptor from optimized complex and ligand from optimized complex. From these energies, the binding free energy is calculated as:

Prime MMGBSA ΔG (Bind) = Complex − Receptor − Ligand

More negative the value of ΔG (Bind), stronger is the binding.

Generating Ligand-Protein Interaction Profile

Ligand protein interaction profile was generated from the pose of ligand showing the highest GlideScore. All the necessary interactions like hydrogen bonding, pi-pi stacking, pi-cation interactions; salt bridges of the HDAC-HDAC inhibitor docked complexes were generated using the default cutoff (4 Å).

Molecular Dynamics Simulation

Molecular dynamics simulation (MDS) is a golden computational technique that helps in monitoring ligand receptor stability and compatibility in an elegant manner. MDS was performed by using Desmond software (Bowers et al., 2006), a relatively novel molecular dynamics (MD) engine by D.E. Shaw Research that can run various molecular simulations including standard MD and simulated annealing (Desmond v4.8) (Shivakumar et al., 2010). The setup files were generated using the system builder option of Desmond. HDAC-HDAC inhibitor complexes were solvated using TIP4P water model (Grover et al., 2012). For specifying shape and size of the repeating unit default orthorhombic boundary conditions were set up. The initial box volume for HDAC2-LAQ824 complex was 463,768 ų and on minimizing, the volume became 416,755 ų. Same parameters for HDAC2-HC-toxin were 468,092 and 419,994 ų respectively. While these values for HDAC1-LAQ824 were 463,596 and 450,033 ų, HC- toxin-HDAC1 displayed these volumes as 463,612 and 450,144 ų. Individual systems were neutralized by adding proper number of counter ions. While HDAC2-LAQ824 and HDAC1-LAQ824 systems were neutralized by adding 2 chloride (Cl⁻) ions and no ions respectively, HDAC2-HC-toxin and HDAC1-HC-toxin were neutralized by adding 2 chloride ions (Cl⁻) and 1 sodium ion (Na⁺). Thus, prior to actual MDS a solvated system including enzyme-ligand complexes as solute and water molecules with oppositely charged ions as solvent was generated. Simulation was carried out under NPT (constant number of atoms, constant pressure and constant temperature) ensemble for 100 ns using the molecular dynamics option of Desmond. Constant temperature (300 K) and pressure (1.01325 bar) was maintained throughout the simulation utilizing Nosé-Hoover thermostat (Hoover, 1985) and Martina-Tobias-Klein method (Martyna et al., 1994; Guo et al., 2010). Moreover, all the solvated systems were relaxed prior to simulation. Detailed information like protein and ligand RMSD, root mean square fluctuation (RMSF) and ligand interaction profile was generated from the simulation trajectory of ligand-receptor complexes using simulation interaction diagram option of Desmond (Pravin et al., 2015).

Creating Hypothesis Using e-Pharmacophores Method

Both structure-based protein-ligand docking and ligand-based pharmacophore modeling are essential parts of drug discovery. While ligand-based technologies are time saving, the structure-based approaches are relatively time consuming but can yield more diverse actives and provide essential target insights. The e-Pharmacophores method incorporates the beneficial aspects of both ligand and structure based approaches by generating structure based pharmacophores which are energetically optimized (Kalyaanamoorthy and Chen, 2013). These pharmacophores can serve as representatives in rapid screening of huge databases. The e-pharmacophores were generated using the auto e-Pharmacophores tool of Phase module of Schrödinger package (Phase v4.9). These pharmacophores were generated by coupling MDS to e-Pharmacophores method via Desmond trajectory clustering. The trajectory from MDS for each HDAC-HDAC inhibitor complex was separated into seven different clusters using the predefined approach. E-pharmacophores were generated from cluster with maximum number of frames which is directly related to ligand stability in active site.

e-Pharmacophores Based Virtual Screening

For validating the e-pharmacophore models, LAQ824 was selected as representative. The e-Pharmacophores of HDAC1-LAQ824 and HDAC2-LAQ824 were used as queries separately and screening was performed against phase database possessing millions of small molecules (Phase v4.9) (Dixon et al., 2006; Natarajan et al., 2016). During virtual screening 50 conformers were generated for each molecule and based on fitness and logic top five hits were selected for each e-Pharmacophores. The hits were docked in extra precision flexible mode (Zhou et al., 2007) against the respective receptors and docking scores were calculated. The docked complexes were subjected to Prime MMGBSA for calculating binding free energy values (Singh and Muthusamy, 2013).

Results and Discussion

RMSD Calculation for Pose Validation

For testing algorithm and selecting docking method pose validation was performed. Self-docking involving the docking of ligands into their native binding site provides a reasonable setup for evaluating docking programs and scoring functions (Sandor et al., 2010). The flexible docking protocol in extra precision mode reproduced a pose very close to crystal pose of ligand. The RMSD value of crystal and redocked pose was found to be 0.39 (Figure 2), authenticating the algorithm for reproducing correct pose (Figure 2). Thus, all the docking studies were performed using extra-precision flexible docking protocol (Jain, 2008).

FIGURE 2

Figure 2. Pose validation for testing docking algorithm. Protein was minimized along with native ligand. Ligand was separated and protein was subjected to grid generation. The separated native ligand was redocked with the host protein using extra precision flexible docking method and root mean square deviation (RMSD) was calculated between native pose (yellow) and redocked pose (dark brown). The native and redocked pose showed RMSD value of 0.39 Å thereby validating the ability of algorithm to determine active site and pose of ligand correctly. Cherry red sphere represents Zinc present at the active site.

Molecular Docking and Evaluation of GScore

Molecular docking, the central tool in drug discovery, provides details regarding the interaction profile of small molecules and a protein at atomic level. Studying such interactions sheds light on the behavior of small molecules in the binding pocket of target proteins. Moreover, these interactions play a considerable role in unraveling the principal biochemical processes (Meng et al., 2011). Accurate structural modeling and correct prediction of activity are regarded as the two chief aims of molecular docking (Ganai, 2016d). In order to gain insights on substantial ligand characteristics arising during the interaction between HDACs and their inhibitors, 5 HDACi from distinct structural groups were chosen (Supplementary Figure 1). These inhibitors with previously determined in vitro IC₅₀ values under a similar experimental setup (Bradner et al., 2010) were docked against HDAC1 and HDAC2 isoforms. The GScore obtained from 10 docked complexes are shown in Figure 3 and Supplementary Table 1. More negative the GScore more is the affinity of ligand toward receptor and vice versa. Among the chosen HDACi, LAQ-824 and other hydroxamates showed most favorable (more negative) GScore. These results are in consistent with the previous reports where hydroxamates and benzamide derivatives showed highest GScore against Class II HDACs and short chain fatty acids like sodium butyrate showed least negative GScore (Ganai et al., 2015b). Similar results were reported by another group against class I HDACs. LAQ824 proved to be the first highest scoring among hydroxamates showing a GScore of −10.74 against HDAC2. While HC-toxin, the cyclic peptide group inhibitor showed GScore of −6.64, valproic acid, the short chain fatty acid group inhibitor showed GScore of −3.61 against HDAC2.

FIGURE 3

Figure 3. GScores of structurally distinct HDACi against HDAC1 and HDAC2 enzymes. More negative the GScore, more is the affinity of ligand toward receptor and vice versa.

However, unlike HDAC2 the first most favorable inhibitor in terms of GScore proved to be LBH-589 (−9.51). This was followed by LAQ824 as evident from a GScore of −9.4. While HC-toxin displayed a GScore of −6.65, valproic acid proved to be feeble inhibitor against HDAC1 as evidenced by least negative GScore (−2.48) respectively. These finding correlate well with the previous findings reporting short chain fatty acids like butyrate and valproic acid as feeble inhibitors of HDACs (Dokmanovic et al., 2007; Rasheed et al., 2007).

Our docking studies revealed that among hydroxamates LAQ824 has strong binding potential especially against the HDAC2 and valproic acid has feeble affinity toward both HDACs. Moreover, HC-toxin, a cyclic peptides showed considerable affinity score against HDAC1 and HDAC2 enzymes. Thus, in the future experiments we tried to explore the reason of differential affinity of hydroxamates, cyclic peptides and short chain fatty acid valproic acid against these Class I isoforms.

Ability to Form Multiple Interactions Markedly Enhances Affinity Score

In order to gain insights about the mechanistic details resulting in the strong affinity score of hydroxamates and the weakest affinity of short chain fatty acid group HDACi, ligand interaction profile of 10 docked complexes were generated. This profile highlights the different interactions possible between ligand and receptor. Hydroxamates like LAQ-824 and LBH589 displayed various interactions with the receptors HDAC2 and HDAC1. These inhibitors only showed highest negative GScore (strong affinity score) during molecular docking studies. Short chain fatty acids, such as phenylbutyrate and valproic acid form only a few interactions with the defined receptors. These inhibitors showed least negative GScore (weak affinity) especially in case of HDAC1. LAQ824, the highest scoring hydroxamate showed three hydrogen bond interactions, one Pi-Pi stacking with HDAC2. While the hydrogen bonding residues with the aforementioned inhibitor were Glu 208, Asp 104 and Gly 154, His 183 was found to be involved in Pi-Pi stacking. On the other hand valproic acid, the short chain fatty acid formed one hydrogen bond with Tyr 308, one salt bridge and metal coordination with active site Zinc ion of HDAC2 (Figure 4A). Similarly hydroxamates LAQ-824 interact with HDAC1 through various bondings in contrast to short chain fatty acid valproate. LAQ-824 interacts via three hydrogen bonds (Gly149, Asp199, and Glu98), four Pi-Pi stacking (Phe205, Phe150, His178, His28) and one Pi-cation interaction (Phe205) with the active site of HDAC1 (Figure 4B). While hydroxamate CRA024781 form salt bridges (Arg39 and Glu103 respectively) with the binding pocket of HDAC2, no such interactions were seen in case of HDAC1 except in case of valproic acid which forms salt bridge with Zinc cation. The only cyclic peptide inhibitor HC-toxin forms two and one hydrogen bonds with the active site residues of HDAC2 (His183 and Tyr308) and HDAC1 (His178) respectively.

FIGURE 4

Figure 4. (A) Ligand interaction diagrams of structurally different HDACi against HDAC2. Diagrams were generated from poses having highest negative GScore using default parameters (upto 4 Å distance from ligand). LAQ824 forms more interactions with the receptor compared to valproic acid and thus shows more negative GScore. Kindly refer lid for more explanation. (B) Ligand interaction profile of different HDACi against HDAC1.

Speaking in general, the high affinity of hydroxamate HDACi may be attributed to higher number of aromatic rings present in these inhibitors. The presence of aromatic rings in these inhibitors results in Pi-Pi stacking and Pi-cation interactions, making their GScore more negative. Moreover, the hydroxamate group in these inhibitors forms hydrogen bonds with HDAC1 and HDAC2 active site besides forming metal coordination with Zinc. Short chain fatty acid valproic acid devoid of any aromatic ring lacks the ability to form Pi-Pi stacking and Pi-cation interactions, which reduces its negative GScore drastically. Taken together, molecular docking studies and ligand interaction profile analysis shows hydroxamates as potent inhibitors against HDAC1 and HDAC2. Moreover, HC-toxin, a cyclic peptide showed better interaction profile against HDAC1 and HDAC2 enzymes. Thus, in nutshell, hydroxamates and cyclic peptide HC-toxin show better interaction profile with these enzymes. In order to confirm these findings we used molecular mechanics generalized born surface area (MMGBSA), an implicit solvation model to calculate the binding free energy of ligand-receptor complexes.

Binding Free Energy Values Aligned Well with in Vitro IC₅₀ Values

The binding free energy of 10 docked complexes was calculated by using MM-GBSA under default conditions (Singh and Muthusamy, 2013). The calculated values showed parallel trend with the in vitro IC₅₀ values of the respective inhibitors (Bradner et al., 2010). More negative the value of binding free energy more is the affinity between ligand and receptor and vice versa. Similarly lower the IC₅₀, more potent is the inhibitor and thus IC₅₀ is inversely proportional to potency of inhibitor (Ganai et al., 2015b). Hydroxamate LAQ824 with an IC₅₀ of 0.003 μM against HDAC2 showed binding free energy (BFE) value of −67.7911 kcal/mol. HC-toxin, a cyclic peptide (IC₅₀ = 0.9 μM) displayed a BFE value of −79.5371 kcal/mol while a value of only −9.09224 kcal/mol was obtained for short chain fatty acid valproic acid (IC₅₀ of 75 μM) against the defined HDAC (Figure 5 and Supplementary Table 2). For HDAC1, LAQ824, a hydroxamate (IC₅₀ 0.001 μM) showed a BFE value of −63.61 kcal/mol while valproic acid with higher value of 51 μM showed lesser negative value of binding free energy (37.09804 kcal/mol) (Bradner et al., 2010). Thus, it is quite evident that hydroxamates having lower IC₅₀ values against HDAC1 and HDAC2 show more negative values of binding free energy while short chain fatty acid valproic acid having higher IC₅₀ value showed very poor value expectedly. The parallel trend between binding free energies and the in vitro IC₅₀ values especially in case of hydroxamates and short chain fatty acid was obtained using the default frozen state where no flexibility was given to receptors HDAC1 or HDAC2. These results are in consistent with the previous findings where maximum correlation has been seen reported for β2-adrenergic receptor agonists and antagonists when the receptor was held frozen (Vilar et al., 2010).

FIGURE 5

Figure 5. BFE values of structurally different HDACi against HDAC1 and HDAC enzymes in frozen state where no flexibility was given to receptor. The binding free energy values reflect the affinity of ligand toward receptor. More negative the value, more is the affinity and vice versa.

Thus calculation of binding free energy using MM-GBSA strongly supported our molecular docking predictions. On analyzing all the three parameters including GScore, interaction profile and binding BFE values we confined our downstream experiments to top scoring hydroxamate LAQ824 and the only cyclic peptide HC-toxin and ignored short chain fatty acid valproic acid due to its poor interaction profile and BFE. The main aim of our study was to generate e-pharmacophores against HDAC1 and HDAC2 enzymes for isoform selective inhibition. Generating these pharmacophores directly from docked complexes does not address the ligand stability and hence the energy issue. Thus, we performed MDS prior to creation of e-Pharmacophores for docked complexes to overcome this impediment.

Molecular Dynamics Simulation Confirmed Receptor Stability throughout Simulation

Though static structure-based approaches like molecular docking and virtual screening have contributed significantly to advanced drug discovery, but they do not take into account, the dynamic nature of proteins (Gunasekaran and Nussinov, 2007). As our studies are related to inhibitor induced conformational changes of active sites where this dynamic nature cannot be breached (Gutteridge and Thornton, 2005; Samsonov et al., 2014), we coupled MDS to these static approaches to overcome the defined loophole. MDS provides in depth insights about small-molecule receptor stability and compatibility and consequently the competency of these molecules to modulate receptor physiological function. Taking these grim facts into consideration, we performed MDS for 100 ns for the most potent inhibitors LAQ824 and HC-toxin in docked state with HDAC1 and HDAC2. RMSD of protein was calculated by aligning all protein frames on the reference frame backbone. Monitoring the RMSD provides crucial details about the structural conformation of protein all through the simulation. Our MDS studies showed protein RMSD well below 3 Å clearly suggesting that docked complexes HDAC1-LAQ824, HDAC2-LAQ824, HDAC1-HC-toxin, and HDAC2-HC-toxin are highly compatible with each other (Figures 6A–D). The RMSF is beneficial for delineating local changes along the protein chain. The secondary structure elements like alpha helices and beta strands portrayed less fluctuation compared to loop regions expectedly (Figures 6A–D). While in case of HDAC1, the residues critical for interaction with LAQ824 are Phe 150, Phe 205, Tyr 303, His 140, Gly 149, Pro 29, and Asp 99 the residues Phe 210, His 145, His 146, Tyr 308 Gly 154, Phe 155, Glu 208, and Asp 104 are crucial in case of HDAC2 (Figures 6E,F). HC-toxin interacts with Tyr 303, Phe 205, Asp 99 and Phe 150 residues of HDAC1 (Figures 6G,H). Regarding HDAC2, HC-toxin portrayed interactions with Tyr 308, Phe 210, Phe 155, Gly 307, and Asp 104 (Figures 7A–H). Previous studies have shown that classical HDACs deacetylate histone substrates using a charge relay mechanism like serine proteases. In this mechanism histidine and aspartate residues participate in addition to Zinc ion (Finnin et al., 1999). Our current findings established that HDACi chelate Zinc ion and thus inactivate classical HDACs (Figure 5). These findings are in consistent with the earlier findings where HDACi have been reported to disrupt charge relay mechanism in case of class II HDACs through Zinc chelation (Ganai et al., 2015b).

FIGURE 6

Figure 6. Molecular dynamics simulation of LAQ824-HDAC1 (A) and LAQ824-HDAC2 (B) docked complexes around 100 ns. RMSD is calculated for all frames and for frame x is:

$RMS{D}_X=\sqrt{\frac{1}{N}{\displaystyle \sum _{i=1}^N{\left(r_i^{\prime }\left(t_x\right))-r_i\left(t_{ref}\right)\right)}^2}}$

Here N represents number of atoms in atom selection; t_ref designates reference time, first frame is selected as reference and is considered as time t = 0; r′ represents the position of selected atoms in frame x after aligning on the reference frame, where frame x is recorded at time t_x. The RMSD of proteins (blue line) is well below 3 clearly suggesting the protein is pretty stable throughout simulation. Moreover, the results show that LAQ824 is quite stable in the active site of these proteins. (C,D) represent the Root Mean Square Fluctuation (RMSF) for LAQ824-HDAC1 and LAQ824-HDAC2. It has importance in characterizing local changes along the protein chain. RMSF for residue i is calculated as:

$RMS{F}_i=\sqrt{\frac{1}{T}{\displaystyle \sum _{i=1}^T<{\left(r_i^{\prime }\left(t\right))-r_i\left(t_{ref}\right)\right)}^2>}}$

In this equation, T signifies trajectory time over which RMSF is calculated, t_ref denotes reference time; r_i represents position of residue i; r′ designates the position of atoms in residue i after superposition on the reference, average of square distance is taken over the selection of atoms in residue is indicated by angle brackets. Peaks represent the areas fluctuating maximum during the simulation. The N and C terminal regions relatively fluctuate more than any other protein part. The loop regions show more fluctuation than rigid protein regions including alpha helices and beta strands. Green colored vertical bars represent the protein residues interacting with ligand (LAQ824). Different types of interactions between LAQ-824-HDAC1 (E,G) and LAQ824-HDAC2 (F,G) throughout the simulation. A value of 1.0 in stacked bar charts represents that during 100% of the simulation time a particular interaction is maintained. Ligand interaction diagrams represent interactions sustaining over 30% of the simulation time. Refer lid given on side for various bonds. It is quite evident that HDACi form metal coordination with Zinc ion present at the active site of HDACs and thus disrupt charge relay mechanism.

FIGURE 7

Figure 7. Diagrams obtained after performing molecular dynamics simulation of HC-toxin-HDAC1 and HC-toxin-HDAC2 docked complexes for around 100 ns. While (A,B) indicate the RMSF of HDAC1 and HDAC2 in docked state with HC-toxin, (C,D) represent the RMSF of the defined HDACs in bound state with the HC-toxin. Different types of interactions between HC-toxin and the two HDACs. (E,G) represent interaction profile of HC-toxin-HDAC1; (F,H) represent same profile for HC-toxin-HDAC2. Consider interactions made by HC-toxin with receptor which includes both HDAC1 and HDAC2 and their cofactor Zn.

Speaking concisely our molecular dynamics studies showed that HDAC1 and HDAC2 are quite stable throughout simulation. Moreover, in the protein-ligand contact diagrams it is quite evident that HDACi target Zinc ion, one of the crucial player in charge relay mechanism. The main aim in this study was to develop e-Pharmacophores from the highly stable pose of ligand which can be obtained only by serious analysis of MDS trajectory. Thus, we performed clustering of MDS output file using Desmond trajectory clustering prior to e-pharmacophore approach.

Creating Pharmacophoric Features of Top Scoring Inhibitors against HDAC1 and HDAC2 Enzymes

Pharmacophore according to Paul Ehrlich is a molecular framework responsible for the biological activity of the drug. According to IUPAC, pharmacophore is an ensemble of steric and electronic features vital for ensuring the optimal supramolecular interactions with a particular biological target and to activate or impede its biological response. The information regarding the common properties among the different binding groups is indispensable for determining the type of inhibitors that are binding the receptor. The pharmacophoric predictions provide such information and thus augment the guidance for the rational design of novel molecules (Zhu et al., 2010). Taking these facts into consideration, we explored the pharmacophoric features of HDACi belonging to structurally distinct groups against HDAC1 and HDAC2 isoforms. In ligand based pharmacophore modeling receptor is not known and the common hypothesis is generated from the active ligands against that receptor (Vuorinen et al., 2014). The energy based models developed in the present work involves the docking of bioactive small-molecule into the active site of receptor. The e-Pharmacophore approach uses the Glide XP scoring function to precisely characterize protein-ligand interactions (Loving et al., 2009; Salam et al., 2009; Ganai, 2016d). As this approach also involves receptor to generate complimentary features, one active compound is sufficient to achieve logical conclusion (Ganai et al., 2015b). As aforementioned, the previous studies have been generated e-pharmacophores directly from docked complexes as thus do not take into account the ligand stability, we used a novel strategy of coupling molecular dynamics to e-pharmacophore via Desmond trajectory clustering, a script which performs hierarchical clustering using MDS output file of protein-ligand complex as input. The MDS trajectory of each docked complex was separated into 7 clusters and the one with maximum number of frames was used for generating e-Pharmacophores. We used auto e-pharmacophore of Phase for generating the critical features of LAQ824, a hydroxamate representative and HC-toxin, a cyclic tetrapeptide against HDAC1 and HDAC2, class I isoforms.

LAQ824 Showed Differential Features against HDAC1 and HDAC2 Enzymes

As aforementioned the pharmacophores were generated from cluster with highest frame frequency. LAQ824 portrayed distinct features against HDAC1 and HDAC2. The defined inhibitor showed five features against HDAC1 namely two hydrogen bond donors (D5) in ZBG and D4 in linker region; two aromatic rings one in linker (R11) and another in cap region (R10); one positive ionizable group (P8) in linker region. Contrary to HDAC1, LAQ824 showed three features against HDAC2 namely one HBD's (D5) in ZBG; two aromatic rings (R9 and R11) in cap and linker region respectively. Thus, LAQ824 showed two features (P8 and D4) lesser against HDAC2 in comparison to HDAC1 (Figures 8A,B).

FIGURE 8

Figure 8. E-Pharmacophoric features of LAQ824 and HC-toxin against HDAC1 and HDAC2 enzymes. LAQ824 shows distinct features against HDAC1 and HDAC2. For HDAC1, LAQ824 portrays five features while for HDAC2 it shows only three features. Two extra features both in linker region of LAQ824 (P8 and D4) were seen for LAQ824 while in the binding pocket of HDAC1 as compared to HDAC2. While HC-toxin portrayed three features (hydrophobic group in linker region, HBA in ZBG and HBD in cap region) against HDAC1, only one feature (HBA in ZBG) was seen against HDAC2. (A) HDAC1- LAQ824, (B) HDAC2- LAQ824, (C) HDAC1- HC-toxin, and (D) HDAC2- HC-toxin.

2 HC-Toxin Showed Two Features Less For Inhibiting HDAC2 Compared to HDAC1

HC-toxin showed three distinct features namely hydrogen bond acceptor (HBA) in ZBG (A1); one hydrophobic group (H13) in linker region and one hydrogen bond donor (HBD) in cap region. For HDAC2, it portrayed only one HBA (A1) in ZBG (Figures 8C,D). Thus, the above findings suggest that for HDAC2 specific inhibition, LAQ824 pharmacophore should have one HBD and one positive ionizable group (both in linker region) less for HDAC2 as compared to HDAC1. While for type-specific inhibition of HDAC1, HC-toxin pharmacophore should have three features located in ZBG, linker and cap, only single feature HBA in ZBG is required for HDAC2 isoform.

LAQ-824 and HC-Toxin Showed Differential High Energy Features against HDAC1 and HDAC2

Energy contribution of each pharmacophoric feature was calculated for LAQ824 and HC-toxin against HDAC1 and HDAC2 enzymes. While positive ionizable group (P8) in linker region (R_L) was maximum scoring (−4.07 kcal/mol) against HDAC1, aromatic ring (R11) in linker (R_L) proved to be best scoring feature (−1.26 kcal/mol) in case of HDAC2. Regarding HC-toxin, HBA in ZBG (A_Z) scored maximum in terms of energy (−0.7 kcal/mol) against HDAC1, followed by hydrophobic group in linker and HBD in cap. However, for HDAC2, HBA in ZBG (A_Z) was the sole top scoring feature (−0.7 kcal/mol). The details of various pharmacophoric features with rank wise scores are summed in Table 1. Summarizing in few words P_L > D_Z > R_L > R_C > D_L in case of LAQ824-HDAC1 while R_L > D_Z > R_C in case of LAQ824-HDAC2 (Table 2). From these studies it is clear that a single inhibitor exhibits differential features against HDAC1 and HDAC2 isoforms and thus in the downstream steps we used the e-Pharmacophores of HDAC1-LAQ824 and HDAC2-LAQ824 as queries in e-Pharmacophore based virtual screening and selected five hits for each of the e-Pharmacophores and tested them against the respective receptors for inhibition.

TABLE 1

Table 1. Ligand interaction profile of HDACi against HDAC1 and HDAC2 enzymes.

TABLE 2

Table 2. Energy contribution of individual pharmacophoric features against HDAC1 and HDAC2, class I isoforms.

Hits Retrieved from e-Pharmacophores Based Virtual Screening Portrayed Promising Affinity toward Therapeutically Relevant HDAC2

The hit molecules retrieved from e-Pharmacophores based virtual screening were selected based on fitness and logic as described in methodology section. For HDAC1, five hits having CACPD2011aCode CACPD2011a-0000302377, CACPD2011a-0001697283, CACPD2011a-0002233975, CACPD2011a-0001697630, and CACPD2011a-0000523061 (Supplementary Figure 2) were finally selected among more than 1000 hits. For HDAC2, CACPD2011a-0001267628, CACPD2011a-0001261600, CACPD2011a-0001267103, CACPD2011a-0000253112, and CACPD2011a-0001261277 (Supplementary Figure 3) were selected among the retrieved hits. The selected hits obtained from HDAC1-LAQ824 and HDAC2-LAQ824 were docked against HDAC1, and HDAC2 respectively. These hits showed negative values of docking score and negative values of binding free energy clearly indicating that these hits do inhibit the respective isoforms and thus validating the e-Pharmacophores hypothesis. Among the selected hits for HDAC1, the first three compounds showed most favorable (more negative) docking as well as BFE values (Figures 9A,B). Regarding the therapeutically relevant HDAC2, all the five hits in general and the hits second, third and fourth in particular portrayed promising affinity toward HDAC2 as evidenced by the expectedly more negative values of docking score and BFE (Figures 9C,D).

FIGURE 9

Figure 9. Docking scores and binding free energy values of virtually screened hits against HDAC1 and HDAC2. While (A,B) shows docking scores and binding free energy values of hits against HDAC1, (C,D) designate the defined parameters for HDAC2.

Conclusion

In this work we have taken the advantage of combinatorial in silico approach including XP-molecular docking, MMGBSA, MDS, trajectory clustering and e-Pharmacophores approach and e-Pharmacophores based virtual screening to exploit the significances of various structural variants in the HDAC inhibitor-HDAC1 and HDAC inhibitor-HDAC2 complexes. HDACi from three distinct structural groups; hydroxamates, cyclic tetrapeptides and short chain fatty acids were docked against HDAC1 and HDAC2 enzymes for evaluating the interaction mechanisms and affinity scores. The predicted BFE values aligned well with the in vitro IC₅₀ values of inhibitors thereby validating the simulation accuracy. Our molecular docking studies demonstrated that hydroxamates like LAQ824, LBH589 etc., and cyclic tetrapeptide HC-toxin display higher affinities with structurally identical HDAC1 and HDAC2, class I isoforms. Moreover, these studies confirmed that inhibitors with higher number of aromatic rings like LAQ824 show enhanced potency while lacking such rings (valproic acid) show feeble affinity toward these enzymes.

Classical HDACs deacetylate histone substrates using charge relay mechanism in which two histidine, two aspartate residues and Zinc ion play a crucial role. Our MDS studies showed that HDACi chelate the Zinc ion by interacting with it thereby disrupting charge relay mechanism, providing further support to the previous findings of Finnin et al. (1999) and Ganai et al. (2015b). We performed MDS of docked complexes followed by MDS trajectory clustering for confirming the stability of protein throughout simulation and for understanding the ligand stability. Our e-Pharmacophores approach proved that HDACi LAQ824 and HC-toxin show differential features against the structurally identical HDAC1 and HDAC2, Class I isoforms. Presence of positive ionizable group and HBD in the linker region of LAQ-824 has significance for inhibiting HDAC1 while absence of these two features is meaningful in case of HDAC2. These findings suggest that presence of positive ionizable group and HBD in linker region of LAQ-824 favors its selectivity toward HDAC1 while absence of these features shifts the selectivity toward HDAC2. The presence of HBD in cap region and hydrophobic group in linker region of HC-toxin promotes its selectivity toward HDAC1 while deficiency of these features enhances its selectivity toward HDAC2.

Quantitative analysis of energy contribution of each pharmacophoric feature was performed for LAQ824-HDAC1, LAQ824-HDAC2, HC-toxin-HDAC1, and HC-toxin-HDAC2 docked complexes. Our studies identified that the positive ionizable group in linker region (P_L) of LAQ824 was the highest scoring pharmacophoric feature against HDAC1 whereas ring in linker (R_L) scores maximum in the binding pocket of HDAC2. Hydrogen bond acceptor (HBA) in ZBG (A_Z) was maximum scoring in case of both HDAC1 and HDAC2. Speaking concisely our in silico screening and e-Pharmacophore models tempt us to speculate that P_L and R_L of LAQ824 have a markedly significant role in inhibiting HDAC1 and HDAC2 enzymes respectively. From the above studies, it is quite evident that e-Pharmacophore models of same inhibitor shows different features against HDAC1 and HDAC2. The hits retrieved from e-Pharmacophores based virtual screening using the HDAC1-LAQ824 and HDAC2-LAQ824 pharmacophores as queries against phase database showed negative values of docking score and BFE thereby validating the accuracy of models. As the studies are focused on HDAC2, hits CACPD2011a-0001261600, CACPD2011a-0001267103, CACPD2011a-0000253112 proved to be most effective in inhibiting HDAC2. These hits with high potency and desired selectivity toward HDAC2 will play a key role in ameliorating neurodegenerative events safely while leaving the function of neuroprotective isoform intact.

Author Contributions

SG and MA designed research and wrote the paper. SG, EA, and RR performed experiments and analyzed data.

Conflict of Interest Statement

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Acknowledgments

The work is supported by Ramanujan Fellowship grant from DST, India and DST-SERB (EMR/2016/000958) to MA. SG acknowledges Financial support from DST-SERB in the form of Start-Up Grant for Young Scientists (YSS/2015/001267), Government of India, and for the procurement of Schrödinger tools (Glide and Prime). RR acknowledges DST INSPIRE, New Delhi, India for doctoral fellowship.

Supplementary Material

The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fnmol.2017.00357/full#supplementary-material

Supplementary Figure 1. Name and structure of HDACi selected for study. The PubChem CID of HDACi is provided in bracket.

Supplementary Figure 2. Structures of five hits obtained through e-Pharmacophores based virtual screening for HDAC1. The 5th structure with CACPD2011aCode (CACPD2011a-0000302377) is the first hit, fourth structure (CACPD2011a-0001697283) is the second hit, third structure (CACPD2011a-0002233975) is the hird hit, structure second is the fourth hit and structure 1 (CACPD20111a-0000523061) is the 5th hit.

Supplementary Figure 3. Structures of hits selected for HDAC2 after virtual screening using LAQ824-HDAC2 e-Pharmacophores as query. Structure second (CACPD2011a-0001267628) is the first hit, structure first is the second hit, structure 5th is the third hit, structures third and fourth are the fourth and 5th hit respectively.

References

Akhtar, M. W., Raingo, J., Nelson, E. D., Montgomery, R. L., Olson, E. N., Kavalali, E. T., et al. (2009). Histone deacetylases 1 and 2 form a developmental switch that controls excitatory synapse maturation and function. J. Neurosci. 29, 8288–8297. doi: 10.1523/JNEUROSCI.0097-09.2009

PubMed Abstract | CrossRef Full Text | Google Scholar

Albani, D., Polito, L., and Forloni, G. (2010). Sirtuins as novel targets for Alzheimer's disease and other neurodegenerative disorders: experimental and genetic evidence. J. Alzheimers Dis. 19, 11–26. doi: 10.3233/JAD-2010-1215

PubMed Abstract | CrossRef Full Text | Google Scholar

Bardai, F. H., Price, V., Zaayman, M., Wang, L., and D'Mello, S. R. (2012). Histone deacetylase-1 (HDAC1) is a molecular switch between neuronal survival and death. J. Biol. Chem. 287, 35444–35453. doi: 10.1074/jbc.M112.394544

PubMed Abstract | CrossRef Full Text | Google Scholar

Bowers, K. J., Chow, E., Xu, H., Dror, R. O., Eastwood, M. P., Gregersen, B. A., et al. (2006). “Scalable algorithms for molecular dynamics simulations on commodity clusters,” in Proceedings of the 2006 ACM/IEEE Conference on Supercomputing, (Tampa, FL: ACM), 84.

Google Scholar

Bradner, J. E., Mak, R., Tanguturi, S. K., Mazitschek, R., Haggarty, S. J., Ross, K., et al. (2010). Chemical genetic strategy identifies histone deacetylase 1 (HDAC1) and HDAC2 as therapeutic targets in sickle cell disease. Proc. Natl. Acad. Sci. U.S.A. 107, 12617–12622. doi: 10.1073/pnas.1006774107

PubMed Abstract | CrossRef Full Text | Google Scholar

Corpet, F. (1988). Multiple sequence alignment with hierarchical clustering. Nucleic Acids Res. 16, 10881–10890. doi: 10.1093/nar/16.22.10881

PubMed Abstract | CrossRef Full Text | Google Scholar

Dasmahapatra, G., Lembersky, D., Kramer, L., Fisher, R. I., Friedberg, J., Dent, P., et al. (2010). The pan-HDAC inhibitor vorinostat potentiates the activity of the proteasome inhibitor carfilzomib in human DLBCL cells in vitro and in vivo. Blood 115, 4478. doi: 10.1182/blood-2009-12-257261

PubMed Abstract | CrossRef Full Text | Google Scholar

Dixon, S. L., Smondyrev, A. M., and Rao, S. N. (2006). PHASE: a novel approach to pharmacophore modeling and 3D database searching. Chem. Biol. Drug Des. 67, 370–372. doi: 10.1111/j.1747-0285.2006.00384.x

PubMed Abstract | CrossRef Full Text | Google Scholar

Dokmanovic, M., Clarke, C., and Marks, P. A. (2007). Histone deacetylase inhibitors: overview and perspectives. Mol. Cancer Res. 5, 981–989. doi: 10.1158/1541-7786.MCR-07-0324

PubMed Abstract | CrossRef Full Text | Google Scholar

Duque-Afonso, J., Yalcin, A., Berg, T., Abdelkarim, M., Heidenreich, O., and Lubbert, M. (2011). The HDAC class I-specific inhibitor entinostat (MS-275) effectively relieves epigenetic silencing of the LAT2 gene mediated by AML1/ETO. Oncogene 30, 3062–3072. doi: 10.1038/onc.2011.32

PubMed Abstract | CrossRef Full Text | Google Scholar

Durham, B. (2012). Novel histone deacetylase (HDAC) inhibitors with improved selectivity for HDAC2 and 3 protect against neural cell death. Biosci. Horiz. 5, 1–7. doi: 10.1093/biohorizons/hzs003

CrossRef Full Text | Google Scholar

d'Ydewalle, C., Bogaert, E., and Van Den Bosch, L. (2012). HDAC6 at the intersection of neuroprotection and neurodegeneration. Traffic 13, 771–779. doi: 10.1111/j.1600-0854.2012.01347.x

PubMed Abstract | CrossRef Full Text | Google Scholar

Eberharter, A., and Becker, P. B. (2002). Histone acetylation: a switch between repressive and permissive chromatin: Second in review series on chromatin dynamics. EMBO Rep. 3, 224–229. doi: 10.1093/embo-reports/kvf053

PubMed Abstract | CrossRef Full Text | Google Scholar

Finnin, M. S., Donigian, J. R., Cohen, A., Richon, V. M., Rifkind, R. A., Marks, P. A., et al. (1999). Structures of a histone deacetylase homologue bound to the TSA and SAHA inhibitors. Nature 401, 188–193.

PubMed Abstract | Google Scholar

Fischle, W., Kiermer, V., Dequiedt, F., and Verdin, E. (2001). The emerging role of class II histone deacetylases. Biochem. Cell Biol. 79, 337–348. doi: 10.1139/o01-116

PubMed Abstract | CrossRef Full Text | Google Scholar

Friesner, R. A., Murphy, R. B., Repasky, M. P., Frye, L. L., Greenwood, J. R., Halgren, T. A., et al. (2006). Extra precision glide: docking and scoring incorporating a model of hydrophobic enclosure for protein-ligand complexes. J. Med. Chem. 49, 6177–6196. doi: 10.1021/jm051256o

PubMed Abstract | CrossRef Full Text | Google Scholar

Ganai, S. A. (2014). “In silico approaches towards safe targeting of class i histone deacetylases,” in Molecular Life Science, eds R. D. Wells, J. S. Bond, J. B. Klinman, S. S. Masters, and E. Bell (New York, NY: Springer Molecular Life Sciences), 1–9.

Google Scholar

Ganai, S. A. (2016a). Histone deacetylase inhibitor pracinostat in doublet therapy: a unique strategy to improve therapeutic efficacy and to tackle herculean cancer chemoresistance. Pharm. Biol. 54, 1926–1935. doi: 10.3109/13880209.2015.1135966

PubMed Abstract | CrossRef Full Text | Google Scholar

Ganai, S. A. (2016b). Histone deacetylase inhibitors modulating non-epigenetic players: the novel molecular targets for therapeutic intervention. Curr. Drug Targets 27, 27. doi: 10.2174/1389450117666160527143257

CrossRef Full Text | Google Scholar

Ganai, S. A. (2016c). Panobinostat: the small molecule metalloenzyme inhibitor with marvelous anticancer activity. Curr. Top. Med. Chem. 16, 427–434. doi: 10.2174/1568026615666150813145800

PubMed Abstract | CrossRef Full Text | Google Scholar

Ganai, S. A. (2016d). Novel approaches towards designing of isoform-selective inhibitors against class ii histone deacetylases: the acute requirement for targetted anticancer therapy. Curr. Top. Med. Chem. 16, 2441–2452. doi: 10.2174/1568026616666160212122609

PubMed Abstract | CrossRef Full Text | Google Scholar

Ganai, S. A. (2017a). Designing isoform-selective inhibitors against Classical HDACs for effective anticancer therapy: insight and perspectives from in silico. Curr. Drug Targets. doi: 10.2174/1389450118666170112130151. [Epub ahead of print].

PubMed Abstract | CrossRef Full Text | Google Scholar

Ganai, S. A. (2017b). Small-molecule modulation of HDAC6 activity: the propitious therapeutic strategy to vanquish neurodegenerative disorders. Curr. Med. Chem. doi: 10.2174/0929867324666170209104030. [Epub ahead of print].

PubMed Abstract | CrossRef Full Text | Google Scholar

Ganai, S. A., Kalladi, S. M., and Mahadevan, V. (2015a). HDAC inhibition through valproic acid modulates the methylation profiles in human embryonic kidney cells. J. Biomol. Struct. Dyn. 33, 1185–1197. doi: 10.1080/07391102.2014.938247

PubMed Abstract | CrossRef Full Text | Google Scholar

Ganai, S. A., Ramadoss, M., and Mahadevan, V. (2016). Histone Deacetylase (HDAC) Inhibitors-emerging roles in neuronal memory, learning, synaptic plasticity and neural regeneration. Curr. Neuropharmacol. 14, 55–71. doi: 10.2174/1570159X13666151021111609

PubMed Abstract | CrossRef Full Text | Google Scholar

Ganai, S. A., Shanmugam, K., and Mahadevan, V. (2015b). Energy-optimised pharmacophore approach to identify potential hotspots during inhibition of Class II HDAC isoforms. J. Biomol. Struct. Dyn. 33, 374–387. doi: 10.1080/07391102.2013.879073

PubMed Abstract | CrossRef Full Text | Google Scholar

Gräff, J., Rei, D., Guan, J.-S., Wang, W.-Y., Seo, J., Hennig, K. M., et al. (2012). An epigenetic blockade of cognitive functions in the neurodegenerating brain. Nature 483, 222–226. doi: 10.1038/nature10849

PubMed Abstract | CrossRef Full Text | Google Scholar

Grover, A., Katiyar, S. P., Jeyakanthan, J., Dubey, V. K., and Sundar, D. (2012). Blocking Protein kinase C signaling pathway: mechanistic insights into the anti-leishmanial activity of prospective herbal drugs from Withania somnifera. BMC Genomics 13:S20. doi: 10.1186/1471-2164-13-S7-S20

PubMed Abstract | CrossRef Full Text | Google Scholar

Guan, J.-S., Haggarty, S. J., Giacometti, E., Dannenberg, J.-H., Joseph, N., Gao, J., et al. (2009). HDAC2 negatively regulates memory formation and synaptic plasticity. Nature 459, 55–60. doi: 10.1038/nature07925

PubMed Abstract | CrossRef Full Text | Google Scholar

Gunasekaran, K., and Nussinov, R. (2007). How different are structurally flexible and rigid binding sites? Sequence and structural features discriminating proteins that do and do not undergo conformational change upon ligand binding. J. Mol. Biol. 365, 257–273. doi: 10.1016/j.jmb.2006.09.062

CrossRef Full Text | Google Scholar

Guo, Z., Mohanty, U., Noehre, J., Sawyer, T. K., Sherman, W., and Krilov, G. (2010). Probing the alpha-helical structural stability of stapled p53 peptides: molecular dynamics simulations and analysis. Chem. Biol. Drug Des. 75, 348–359. doi: 10.1111/j.1747-0285.2010.00951.x

PubMed Abstract | CrossRef Full Text | Google Scholar

Gutteridge, A., and Thornton, J. (2005). Conformational changes observed in enzyme crystal structures upon substrate binding. J. Mol. Biol. 346, 21–28. doi: 10.1016/j.jmb.2004.11.013

PubMed Abstract | CrossRef Full Text | Google Scholar

Hoover, W. G. (1985). Canonical dynamics: equilibrium phase-space distributions. Phys. Rev. A Gen. Phys. 31, 1695–1697. doi: 10.1103/PhysRevA.31.1695

PubMed Abstract | CrossRef Full Text | Google Scholar

Jain, A. N. (2008). Bias, reporting, and sharing: computational evaluations of docking methods. J. Comput. Aided Mol. Des. 22, 201–212. doi: 10.1007/s10822-007-9151-x

PubMed Abstract | CrossRef Full Text | Google Scholar

Kalyaanamoorthy, S., and Chen, Y. P. (2013). Energy based pharmacophore mapping of HDAC inhibitors against class I HDAC enzymes. Biochim. Biophys. Acta 1834, 317–328. doi: 10.1016/j.bbapap.2012.08.009

PubMed Abstract | CrossRef Full Text | Google Scholar

Lauffer, B. E., Mintzer, R., Fong, R., Mukund, S., Tam, C., Zilberleyb, I., et al. (2013). Histone deacetylase (HDAC) inhibitor kinetic rate constants correlate with cellular histone acetylation but not transcription and cell viability. J. Biol. Chem. 288, 26926–26943. doi: 10.1074/jbc.M113.490706

CrossRef Full Text | Google Scholar

Lee, J. H., Yao, Y., Mahendran, A., Ngo, L., Venta-Perez, G., Choy, M. L., et al. (2015). Creation of a histone deacetylase 6 inhibitor and its biological effects [corrected]. Proc. Natl. Acad. Sci. U.S.A. 112, 12005–12010. doi: 10.1073/pnas.1515882112

PubMed Abstract | CrossRef Full Text | Google Scholar

Loving, K., Salam, N. K., and Sherman, W. (2009). Energetic analysis of fragment docking and application to structure-based pharmacophore hypothesis generation. J. Comput. Aided Mol. Des. 23, 541–554. doi: 10.1007/s10822-009-9268-1

PubMed Abstract | CrossRef Full Text | Google Scholar

Lyne, P. D., Lamb, M. L., and Saeh, J. C. (2006). Accurate prediction of the relative potencies of members of a series of kinase inhibitors using molecular docking and MM-GBSA scoring. J. Med. Chem. 49, 4805–4808. doi: 10.1021/jm060522a

PubMed Abstract | CrossRef Full Text | Google Scholar

Martyna, G. J., Tobias, D. J., and Klein, M. L. (1994). Constant pressure molecular dynamics algorithms. J. Chem. Phys. 101, 4177–4189. doi: 10.1063/1.467468

CrossRef Full Text | Google Scholar

Meng, X.-Y., Zhang, H.-X., Mezei, M., and Cui, M. (2011). Molecular docking: a powerful approach for structure-based drug discovery. Curr. Comput. Aided Drug Des. 7, 146–157. doi: 10.2174/157340911795677602

PubMed Abstract | CrossRef Full Text | Google Scholar

Millard, C. J., Watson, P. J., Celardo, I., Gordiyenko, Y., Cowley, S. M., Robinson, C. V., et al. (2013). Class I HDACs share a common mechanism of regulation by inositol phosphates. Mol. Cell 51, 57–67. doi: 10.1016/j.molcel.2013.05.020

PubMed Abstract | CrossRef Full Text | Google Scholar

Morris, M. J., Karra, A. S., and Monteggia, L. M. (2010). Histone deacetylases govern cellular mechanisms underlying behavioral and synaptic plasticity in the developing and adult brain. Behav. Pharmacol. 21, 409–419. doi: 10.1097/FBP.0b013e32833c20c0

PubMed Abstract | CrossRef Full Text | Google Scholar

Morris, M. J., Mahgoub, M., Na, E. S., Pranav, H., and Monteggia, L. M. (2013). Loss of histone deacetylase 2 improves working memory and accelerates extinction learning. J. Neurosci. 33, 6401–6411. doi: 10.1523/JNEUROSCI.1001-12.2013

PubMed Abstract | CrossRef Full Text | Google Scholar

Mottamal, M., Zheng, S., Huang, T. L., and Wang, G. (2015). Histone deacetylase inhibitors in clinical studies as templates for new anticancer agents. Molecules 20, 3898–3941. doi: 10.3390/molecules20033898

PubMed Abstract | CrossRef Full Text | Google Scholar

Natarajan, P., Priyadarshini, V., Pradhan, D., Manne, M., Swargam, S., Kanipakam, H., et al. (2016). E-pharmacophore-based virtual screening to identify GSK-3β inhibitors. J. Recept. Signal Transduct. Res. 36, 445–458. doi: 10.3109/10799893.2015.1122043

PubMed Abstract | CrossRef Full Text | Google Scholar

Pravin, N., Devaraji, V., and Raman, N. (2015). Targeting protein kinase and DNA molecules by diimine-phthalate complexes in antiproliferative activity. Int. J. Biol. Macromol. 79, 837–855. doi: 10.1016/j.ijbiomac.2015.06.001

PubMed Abstract | CrossRef Full Text | Google Scholar

Qi, H. Y., Zhang, Z. J., Li, Y. J., and Fang, X. D. (2011). [Role of chromatin conformation in eukaryotic gene regulation]. Yi Chuan 33, 1291–1299. doi: 10.3724/SP.J.1005.2011.01291

PubMed Abstract | CrossRef Full Text | Google Scholar

Rasheed, W. K., Johnstone, R. W., and Prince, H. M. (2007). Histone deacetylase inhibitors in cancer therapy. Expert Opin. Investig. Drugs 16, 659–678. doi: 10.1517/13543784.16.5.659

PubMed Abstract | CrossRef Full Text | Google Scholar

Rivieccio, M. A., Brochier, C., Willis, D. E., Walker, B. A., D'Annibale, M. A., McLaughlin, K., et al. (2009). HDAC6 is a target for protection and regeneration following injury in the nervous system. Proc. Natl. Acad. Sci. U.S.A. 106, 19599–19604. doi: 10.1073/pnas.0907935106

PubMed Abstract | CrossRef Full Text | Google Scholar

Ropero, S., and Esteller, M. (2007). The role of histone deacetylases (HDACs) in human cancer. Mol. Oncol. 1, 19–25. doi: 10.1016/j.molonc.2007.01.001

PubMed Abstract | CrossRef Full Text | Google Scholar

Saha, R. N., and Pahan, K. (2006). HATs and HDACs in neurodegeneration: a tale of disconcerted acetylation homeostasis. Cell Death Differ. 13, 539–550. doi: 10.1038/sj.cdd.4401769

PubMed Abstract | CrossRef Full Text | Google Scholar

Saiyed, Z. M., Gandhi, N., Agudelo, M., Napuri, J., Samikkannu, T., Reddy, P. V., et al. (2011). HIV-1 Tat upregulates expression of histone deacetylase-2 (HDAC2) in human neurons: implication for HIV-associated neurocognitive disorder (HAND). Neurochem. Int. 58, 656–664. doi: 10.1016/j.neuint.2011.02.004

PubMed Abstract | CrossRef Full Text | Google Scholar

Salam, N. K., Nuti, R., and Sherman, W. (2009). Novel method for generating structure-based pharmacophores using energetic analysis. J. Chem. Inf. Model. 49, 2356–2368. doi: 10.1021/ci900212v

PubMed Abstract | CrossRef Full Text | Google Scholar

Samsonov, S. A., Gehrcke, J. P., and Pisabarro, M. T. (2014). Flexibility and explicit solvent in molecular-dynamics-based docking of protein-glycosaminoglycan systems. J. Chem. Inf. Model. 54, 582–592. doi: 10.1021/ci4006047

PubMed Abstract | CrossRef Full Text | Google Scholar

Sandor, M., Kiss, R., and Keseru, G. M. (2010). Virtual fragment docking by Glide: a validation study on 190 protein-fragment complexes. J. Chem. Inf. Model. 50, 1165–1172. doi: 10.1021/ci1000407

PubMed Abstract | CrossRef Full Text | Google Scholar

Sastry, G. M., Adzhigirey, M., Day, T., Annabhimoju, R., and Sherman, W. (2013). Protein and ligand preparation: parameters, protocols, and influence on virtual screening enrichments. J. Comput. Aided Mol. Des. 27, 221–234. doi: 10.1007/s10822-013-9644-8

PubMed Abstract | CrossRef Full Text | Google Scholar

Seto, E., and Yoshida, M. (2014). Erasers of histone acetylation: the histone deacetylase enzymes. Cold Spring Harb. Perspect. Biol. 6:a018713. doi: 10.1101/cshperspect.a018713

PubMed Abstract | CrossRef Full Text | Google Scholar

Shivakumar, D., Williams, J., Wu, Y., Damm, W., Shelley, J., and Sherman, W. (2010). Prediction of absolute solvation free energies using molecular dynamics free energy perturbation and the opls force field. J. Chem. Theory Comput. 6, 1509–1519. doi: 10.1021/ct900587b

PubMed Abstract | CrossRef Full Text | Google Scholar

Singh, K. D., and Muthusamy, K. (2013). Molecular modeling, quantum polarized ligand docking and structure-based 3D-QSAR analysis of the imidazole series as dual AT1 and ETA receptor antagonists. Acta Pharmacol. Sin. 34, 1592–1606. doi: 10.1038/aps.2013.129

CrossRef Full Text | Google Scholar

Sinha, S., Goyal, S., Somvanshi, P., and Grover, A. (2016). Mechanistic insights into the binding of Class IIa HDAC inhibitors toward spinocerebellar ataxia type-2: a 3D-QSAR and pharmacophore modeling approach. Front. Neurosci. 10:606. doi: 10.3389/fnins.2016.00606

PubMed Abstract | CrossRef Full Text | Google Scholar

Steinbrecher, T., Zhu, C., Wang, L., Abel, R., Negron, C., Pearlman, D., et al. (2017). Predicting the effect of amino acid single-point mutations on protein stability—large-scale validation of MD-based relative free energy calculations. J. Mol. Biol. 429, 948–963. doi: 10.1016/j.jmb.2016.12.007

PubMed Abstract | CrossRef Full Text | Google Scholar

Van Den Driessche, G., and Fourches, D. (2017). Adverse drug reactions triggered by the common HLA-B*57:01 variant: a molecular docking study. J. Cheminform. 9:13. doi: 10.1186/s13321-017-0202-6

PubMed Abstract | CrossRef Full Text | Google Scholar

Ververis, K., Hiong, A., Karagiannis, T. C., and Licciardi, P. V. (2013). Histone deacetylase inhibitors (HDACIs): multitargeted anticancer agents. Biologics 7, 47–60. doi: 10.2147/BTT.S29965

PubMed Abstract | CrossRef Full Text | Google Scholar

Vilar, S., Karpiak, J., and Costanzi, S. (2010). Ligand and structure-based models for the prediction of ligand-receptor affinities and virtual screenings: development and application to the β₂-adrenergic receptor. J. Comput. Chem. 31, 707–720. doi: 10.1002/jcc.21346

CrossRef Full Text | Google Scholar

Volmar, C.-H., and Wahlestedt, C. (2015). Histone deacetylases (HDACs) and brain function. Neuroepigenetics 1, 20–27. doi: 10.1016/j.nepig.2014.10.002

CrossRef Full Text | Google Scholar

Vuorinen, A., Engeli, R., Meyer, A., Bachmann, F., Griesser, U. J., Schuster, D., et al. (2014). Ligand-based pharmacophore modeling and virtual screening for the discovery of novel 17β-hydroxysteroid dehydrogenase 2 inhibitors. J. Med. Chem. 57, 5995–6007. doi: 10.1021/jm5004914

PubMed Abstract | CrossRef Full Text | Google Scholar

Wagner, F. F., Zhang, Y. L., Fass, D. M., Joseph, N., Gale, J. P., Weïwer, M., et al. (2015). Kinetically selective inhibitors of Histone Deacetylase 2 (HDAC2) as cognition enhancers. Chem. Sci. 6, 804–815. doi: 10.1039/C4SC02130D

PubMed Abstract | CrossRef Full Text | Google Scholar

Watroba, M., and Szukiewicz, D. (2016). The role of sirtuins in aging and age-related diseases. Adv. Med. Sci. 61, 52–62. doi: 10.1016/j.advms.2015.09.003

PubMed Abstract | CrossRef Full Text | Google Scholar

Welberg, L. (2012). Neurodegeneration: HDAC2 mediates cognitive impairment. Nat. Rev. Neurosci. 13, 223–223. doi: 10.1038/nrn3222

CrossRef Full Text | Google Scholar

Witt, O., Deubzer, H. E., Milde, T., and Oehme, I. (2009). HDAC family: what are the cancer relevant targets? Cancer Lett. 277, 8–21. doi: 10.1016/j.canlet.2008.08.016

PubMed Abstract | CrossRef Full Text | Google Scholar

Yang, X. J., and Seto, E. (2007). HATs and HDACs: from structure, function and regulation to novel strategies for therapy and prevention. Oncogene 26, 5310–5318. doi: 10.1038/sj.onc.1210599

PubMed Abstract | CrossRef Full Text | Google Scholar

Yang, X. J., and Seto, E. (2008). Lysine acetylation: codified crosstalk with other posttranslational modifications. Mol. Cell 31, 449–461. doi: 10.1016/j.molcel.2008.07.002

PubMed Abstract | CrossRef Full Text | Google Scholar

Zhou, J., Li, M., Chen, N., Wang, S., Luo, H.-B., Zhang, Y., et al. (2015). Computational design of a time-dependent histone deacetylase 2 selective inhibitor. ACS Chem. Biol. 10, 687–692. doi: 10.1021/cb500767c

PubMed Abstract | CrossRef Full Text | Google Scholar

Zhou, Z., Felts, A. K., Friesner, R. A., and Levy, R. M. (2007). Comparative performance of several flexible docking programs and scoring functions: enrichment studies for a diverse set of pharmaceutically relevant targets. J. Chem. Inf. Model. 47, 1599–1608. doi: 10.1021/ci7000346

PubMed Abstract | CrossRef Full Text | Google Scholar

Zhu, Y., Li, H. F., Lu, S., Zheng, Y. X., Wu, Z., Tang, W. F., et al. (2010). Investigation on the isoform selectivity of histone deacetylase inhibitors using chemical feature based pharmacophore and docking approaches. Eur. J. Med. Chem. 45, 1777–1791. doi: 10.1016/j.ejmech.2010.01.010

PubMed Abstract | CrossRef Full Text | Google Scholar