Docking and binding free energy calculations of sirtuin inhibitors

Docking and binding free energy calculations of sirtuin inhibitors

Accepted Manuscript Docking and binding free energy calculations of sirtuin inhibitors Berin Karaman, Wolfgang Sippl PII: S0223-5234(15)00141-5 DOI:...

4MB Sizes 0 Downloads 77 Views

Accepted Manuscript Docking and binding free energy calculations of sirtuin inhibitors Berin Karaman, Wolfgang Sippl PII:

S0223-5234(15)00141-5

DOI:

10.1016/j.ejmech.2015.02.045

Reference:

EJMECH 7727

To appear in:

European Journal of Medicinal Chemistry

Received Date: 27 November 2014 Revised Date:

25 January 2015

Accepted Date: 22 February 2015

Please cite this article as: B. Karaman, W. Sippl, Docking and binding free energy calculations of sirtuin inhibitors, European Journal of Medicinal Chemistry (2015), doi: 10.1016/j.ejmech.2015.02.045. This is a PDF file of an unedited manuscript that has been accepted for publication. As a service to our customers we are providing this early version of the manuscript. The manuscript will undergo copyediting, typesetting, and review of the resulting proof before it is published in its final form. Please note that during the production process errors may be discovered which could affect the content, and all legal disclaimers that apply to the journal pertain.

ACCEPTED MANUSCRIPT

AC C

EP

TE D

M AN U

SC

RI PT

GRAPHICAL ABSTRACT

ACCEPTED MANUSCRIPT

Docking and binding free energy calculations of

Institute of Pharmacy, Martin-Luther-University of Halle-Wittenberg, 06120 Halle

M AN U



SC

Berin Karaman† and Wolfgang Sippl†

RI PT

sirtuin inhibitors

AC C

EP

TE D

(Saale), Germany

ACCEPTED MANUSCRIPT ABSTRACT Sirtuins form a unique and highly conserved class of NAD+-dependent lysine deacylases. Among these the human subtypes Sirt1-3 has been implicated in the pathogenesis of numerous diseases such as cancer, metabolic syndromes, viral diseases and neurological disorders. Most of the sirtuin

RI PT

inhibitors that have been identified so far show limited potency and/or isoform selectivity. Here, we introduce a promising method to generate protein-inhibitor complexes of human Sirt1, Sirt2 and Sirt3 by means of ligand docking and molecular dynamics simulations. This method highly reduces

SC

the complexity of such applications and can be applied to other protein targets beside sirtuins. To the best of our knowledge, we present the first binding free energy method developed by using a

M AN U

validated data set of sirtuin inhibitors, where both a fair number of compounds (33 thieno[3,2d]pyrimidine-6-carboxamide derivatives) was developed and tested in the same laboratory and also crystal structures in complex with the enzyme have been reported. A significant correlation between binding free energies derived from MM-GBSA calculations and in vitro data was found for all three sirtuin subtypes. The developed MM-GBSA protocol is computationally inexpensive and can be

TE D

applied as a post-docking filter in virtual screening to find novel Sirt1-3 inhibitors as well as to

KEYWORDS

EP

prioritize compounds with similar chemical structures for further biological characterization.

AC C

Sirtuin, MM-GBSA, MM-PBSA, virtual screening.

ACCEPTED MANUSCRIPT 1. INTRODUCTION Reliable estimation of protein-inhibitor binding free energies represents a holy grail in drug design. Over the last few years, different approaches for predicting the biological activity of small mole-

RI PT

cules, each with individual trade-offs of accuracy versus computation time, have been developed [1-8]. Often, X-ray structures are not available for all ligands under study. Therefore, molecular docking is applied to predict the binding pose of the ligands. The two major tasks of molecular docking include the determination of the correct orientation within the binding site and the accurate

SC

assessment of the protein-ligand affinity. Whereas docking methods already show satisfying results for correctly predicting of the binding modes, the second task still remains a key challenge in com-

M AN U

putational chemistry [2,9]. Considering the balance between the accuracy of calculations and the computational cost, scoring functions were developed on the basis of simplified empirical force fields or potentials of mean force using certain approximations, e.g. they usually do not take solvation effects or protein flexibility into account [4,5,10]. Thus, accurate prediction of the binding af-

methods for its evaluation.

TE D

finity using docking methodologies remains an elusive goal and there is the need of more precise

Correct ranking of compounds according to their binding affinity is another critical issue. Beside

EP

traditional scoring methods, more sophisticated and demanding calculations of binding free energy (BFE) such as linear interaction energy (LIE) [11], molecular mechanic/Poisson-Boltzmann (Gen-

AC C

eralized Born) surface area (MM-PB(GB)SA) [12], free-energy perturbation (FEP) and thermodynamic integration (TI) [13,14] have been developed. Although some studies show successful examples of lead optimization using FEP and TI methods [15-17], they are still rarely applied in the drug discovery process, mainly due to the high computational costs. Indeed, the success of exhilarating the use of these techniques (especially in industrial settings) depends mainly on the minimization of human resources for data analysis. Nevertheless, the recent increase in computing power made it possible to apply BFE calculation methods on quite large compound datasets. This allows using them not only in the late stages of drug discovery, but also as a tool for post-processing of vir-

ACCEPTED MANUSCRIPT tual screening results. The number of studies, reporting on the successful application of MM-PBSA, MM-GBSA methods for the estimation of protein-ligand binding affinities, is increasing, and automated procedures for BFE estimation applicable for a large compound selection have been developed [18-24]. Generally, satisfactory correlations have been obtained between the calculated bind-

RI PT

ing affinities and the experimental data especially for ligand series with high structural similarity [17,20,25]. Although, previous applications of MM-PBSA and MM-GBSA proved to be very useful and efficient in numerous macromolecular analyses, there remains still room for improvement in

SC

many aspects of existing methods, including the estimation of the solute entropy and the solvation free energy of polar compounds and of charged groups, the sampling of conformational space and

M AN U

the parameters used for the description of molecular systems [3,17].

NAD+-dependent protein deacylases (sirtuins) are a unique class of enzymes that are conserved from bacteria to humans. Initially recognized as protein deacetylases, they have lately been shown to catalyze further posttranslational modifications such as demyristoylation [26], demalonylation, desuccinylation [27] or ADP-ribosylation [28]. Through a multitude of protein substrates they are

TE D

involved in key cellular processes including metabolic sensing, regulation of mitosis and aging. Seven sirtuin histone deacetylase subtypes (Sirt1-7) have been identified in mammals which share highly conserved catalytic NAD+ and acetyl-lysine binding sites [29,30]. The human sirtuin sub-

EP

types Sirt1-3 represent interesting targets for the treatment of age related diseases, neurological disorders (like Parkinson`s and Alzheimer`s diseases), metabolic syndromes (such as diabetes and obe-

AC C

sity), viral diseases and cancer [31-40]. The biochemistry of sirtuins has been studied with great effort and three-dimensional structures of several human isoforms were solved that give insight into the different stages of the catalytic cycle [41-52]. Nevertheless, a structure-derived rationale for an isoform-selective inhibition remained elusive. The highly conserved amino acid sequence and the high structural similarity of the catalytic core between the members of the sirtuin family makes it challenging to come up with selective drug-like inhibitors. Nevertheless, we and others have reported several series of sirtuin inhibitors that might be used as starting point to design more selective

ACCEPTED MANUSCRIPT compounds [32,53-55]. Due to the high flexibility of the substrate and cofactor binding sites observed in sirtuins, the structure based optimization of inhibitors is still a problem. Recently, Disch et al. reported a series of thieno[3,2-d]pyrimidine-6-carboxamide derivatives, which showed inhibition of Sirt1-3 [43]. For this series of compounds crystal structures in complex with Sirt3 could be re-

ploration of structure activity relationships.

RI PT

solved that revealed the binding mode of this class of compounds and also enabled the further ex-

In this study, we tested whether the MM-GBSA approach is able to provide predictive models for

SC

the series of published Sirt1-3 inhibitors [43]. The obtained results show that binding free energy calculations using a small number of snapshots retrieved from MD simulations give a significant

M AN U

correlation with the experimental data. The established protocol appears to be useful for discovering novel sirtuin inhibitors and can be applied as post-docking filter for virtual screenings of larger compound databases.

1.1. Inhibitor Data Set

TE D

MATERIALS AND METHODS

A data set of 33 thieno[3,2-d]pyrimidine-6-carboxamide derivatives tested against human Sirt1-3

EP

was used in the current study [43]. Structural information was downloaded from EMBL-EBI (European Molecular Biology Laboratory- European Bioinformatics Institute) website [56].

AC C

Canonical smiles were used to build up the 3D structures of the molecules, and a subsequent energy minimization was applied using MMFF94x force-field implemented in Molecular Operating Environment System (MOE) 2012.10 [57]. The chemical structures of the studied inhibitors are displayed in Chart 1. Compounds are named using the same numbering as reported by Disch et al. [43].

1.2. Preparation of Protein-Inhibitor Complexes 1.2.1. Protein-Inhibitor Complexes for Sirt3

ACCEPTED MANUSCRIPT The reported compounds were shown to inhibit Sirt3 with an IC50 in the range of >50 - 0.0032 µM. X-ray crystal structures of Sirt3 in complex with the three inhibitors 11c, 28 and 31 (Chart 1) were downloaded from the Protein Data Bank (PDB) [58,59] (PDB ID; 4JSR, 4JT8 and 4JT9, respectively) [43]. All protein structures were prepared by using the Structure Preparation module

RI PT

in MOE 2012.10 [57]. Protonation states of the structures were determined using the Protonate 3D module in MOE. Protein structures were energy minimized using the AMBER99 force field [60]

SC

using a tethering force constant of (3/2) kT / σ2 (σ = 0.5 Å) for all atoms during the minimization.

1.2.2. Protein-Inhibitor Complexes Generated for Sirt1 and Sirt2

M AN U

The 33 thieno[3,2-d]pyrimidine-6-carboxamide derivatives were also tested in vitro on Sirt1 and Sirt2 sirtuin subtypes and were shown to inhibit Sirt1 (with an IC50 in the range of 50 - 0.0036 µM) and Sirt2 (with an IC50 in the range of 50 - 0.0018 µM) with similar range of IC50 values obtained for Sirt3 inhibition [43]. The activity data of 28 compounds against Sirt1 and Sirt3 and 30

TE D

compounds against Sirt2 were reported [43]. However, no inhibitor bound X-ray structure was resolved for Sirt1 and Sirt2. Yet, pocket residues (within 4.5 Å of superposed inhibitor 31, found in the Sirt3-31 crystal structure) of Sirt3 share a 100 % sequence identity with the isoforms Sirt1 and

EP

Sirt2. In addition, the inhibitors share the same mechanism of inhibition [61]. These facts point out that the compounds most likely bind Sirt1 and Sirt2 in the same manner seen in Sirt3 inhibitor

AC C

bound crystal structure. Available crystal structures of human Sirt1 (PDB ID 4IG9; Sirt1-apo, 4KXQ; Sirt1-ADPR [49] and 4I5I; Sirt1-NAD+-Ex527 analog [48]) and Sirt2 (PDB ID 3ZGO; Sirt2-apo and 3ZGV; Sirt2-ADPR) [45] were taken from the PDB. Protein structures were prepared as mentioned before. Sirtuin structures were superimposed on their backbone atoms using the Superpose module in MOE and using the pocket residues surrounding the inhibitor found in the Sirt3-31 crystal structure. The closest conformation of Sirt1 and Sirt2 structures to the inhibitor bound Sirt3 structures were chosen by taking into account the Root Mean Square Deviation (RMSD) values (Supporting Information Figure S2 and S4) and a careful visual analysis of the

ACCEPTED MANUSCRIPT pocket residues` placement and orientation. We carried out molecular docking studies of the 33 thieno[3,2-d]pyrimidine-6-carboxamide derivatives on Sirt1 and Sirt2 to test whether the selected structures were able to get a binding mode for the inhibitors as seen in the Sirt3 crystal structures. In addition, we applied a short (1ns) molecular dynamics (MD) simulation using AMBER 12 [62]

RI PT

to relax the Sirt1 and Sirt2 structures in complex with inhibitor. Finally, the stability of the proteininhibitor complexes, generated for Sirt1 and Sirt2, were examined by means of a longer (20ns) MD

SC

simulation.

M AN U

1.3. Molecular Docking

Protein structures were prepared as described before. All molecules except the zinc ion were removed from the structures before docking. Structural bridging water molecules were included in the binding site of the protein structures before docking where it is mentioned. Docking studies were performed using the GLIDE program (Schrödinger Suite 2012-5.8) [63-66]. The position of

TE D

the inhibitors in the crystal structures was used to define the grid box location (10 Å radius) and then the box length was increased to 14 Å radius in the Y direction. The dockings were performed using the GLIDE extra precision (GLIDE-XP) mode treating the ligand as flexible. 20 docking

EP

poses were calculated for each inhibitor. All other options were left at their default values. The topranked docking poses were visually analyzed together with the protein structure using the program

AC C

MOE.

1.4. Molecular Dynamics (MD) Simulation MD simulations were carried out using AMBER 12 [62] and the AMBER 2003 force field [67,68]. Atom types and AM1-BCC atomic charges [69,70] were generated for each ligand using the Antechamber module. Ligand parameters were taken from the general AMBER force field GAFF [71]. Preparation of the ligand-protein complex, addition of counter ions, solvation, and preparation

ACCEPTED MANUSCRIPT of parameter/topology and coordinate files were carried out using the LEaP module in AMBER. Parameters and libraries for zinc binding residues were defined as described by Pang Lab [72]. The system was solvated using the water model TIP3BOX [73] and a margin of 10 Å. Two consecutive steps of minimization were carried out. In the first step 3000 iterations (first 1000 steepest descent

RI PT

and then 2000 conjugate gradient) and in the second step 4000 iterations (first 2000 steepest descent and then 2000 conjugate gradient) were applied to the system. In the first step, atom coordinates for the amino acid residues, ligand atoms and conserved water molecules were restrained to their initial

SC

coordinates with a force constant of 500 kcal mol-1Å-2 to relieve the bad van der Waals (vdW) contacts in the surrounding solvent and thus to minimize the positions of the inserted water

M AN U

molecules and ions. In the second step, restraints on atoms were removed and the whole system was minimized freely to relieve bad contacts in the entire system. The system was then equilibrated at 300 K through 100 ps of MD with a time step of 2 fs per step. A constant volume periodic boundary was set to equilibrate the temperature of the system by the Langevin dynamics [74] using a collision frequency of 1 ps−1 during the temperature equilibration routine. The protein and ligand atoms were

TE D

restrained with a weak force constant of 10 kcal mol-1Å-2. The final coordinates obtained after the temperature equilibration step were then used for a longer MD routine either 1ns or 20ns MD simulation is applied where it is stated) during which the temperature was kept at 300 K by the

EP

Langevin dynamics using a collision frequency of 1 ps-1. Constant pressure periodic boundary was used to maintain the pressure of the system at 1 bar using isotropic pressure scaling with a

AC C

relaxation time of 2 ps. During the temperature equilibration and MD routines a non-bonded cut-off distance of 10.0 Å was used applying the Particle Mesh Ewald (PME) method [75] for calculating the full electrostatic energy of the periodic system and the SHAKE algorithm [76] was adjusted to constrain all bonds involving hydrogen.

1.5. Binding Free Energy Calculations Binding free energies for the inhibitors under study were calculated using the top-ranked docking

ACCEPTED MANUSCRIPT poses obtained for Sirt1, 2 and 3. Structurally conserved water molecules included for the docking studies were maintained during the geometry optimization of the complexes, but were removed for the binding free energy calculations. The protein-inhibitor complexes were first energy minimized and then equilibrated using a 100 ps MD run and the same parameters as described before. Three

RI PT

different systems were considered for the binding free energy calculations: a) the protein-inhibitor complex which was energy minimized in explicit water, b) 10 snapshots for each complex extracted evenly from the MD simulation and c) the final complex from the MD equilibration phase. The

SC

binding free energy was then estimated using either the MM-PBSA or MM-GBSA module implemented in AMBER 12.

M AN U

The MM-PB(GB)SA approach combines molecular mechanics and continuum solvent models to predict the protein-ligand binding free energy. The binding strength of a ligand can be calculated as the difference between the free energy of the complex and the sum of the free energies of its components; the receptor and ligand (eq1).

TE D

(eq1) ∆Gbind = Gcomplex – (Greceptor + Gligand) The total free energy of each component in eq1 can be computed as the sum of the gas-phase energy (EMM), the solvation free energy (∆Gsol) and the configurational entropy contributions (-T∆S)

EP

(eq2).

(eq2) Gmolecule = (EMM ) + ( ∆Gsol ) – (T∆S)

AC C

EMM energy is the molecular mechanics energy of the molecule that includes Eint (internal), Eele (electrostatic), Evdw (van der Waals) energies. The first EMM energy term; Eint treats bond, angle and dihedral energies, whereas the second EMM energy term Eele accounts for the Coulombic interactions between atoms occurring as a result of positive and negative partial atomic charges and the Evdw term takes into account short and long ranged interactions among atoms. The component energy terms of EMM energy were calculated using the sander module of AMBER. Gsol is the solvation energy term expressed as the sum of electrostatic solvation energy (polar contribution), GPB/GB, and the non-electrostatic solvation energy (nonpolar contribution), GSA (eq3). The polar solvation

ACCEPTED MANUSCRIPT energy contribution was calculated by solving either Poisson-Boltzmann (PB) equation (MM-PBSA calculation) with the PBSA module or its variant Generalized Born (GB) equation (MM-GBSA calculation) with the MM-GBSA module in AMBER. Different GB models have been developed and implemented in AMBER such as GBHCT (igb=1), GBOC1 (igb=2), GBOC2 (igb=5), and the more

RI PT

recently derived GBn models (igb=7, 8) [62,77-81]. In the current work we compared different MM-GBSA solvation models as well as the MM-PBSA method. The value of the implicit solvent dielectric constant and the solute dielectric constant for GB calculations were set to 80 and 1,

SC

respectively. The solvent probe radius was set to 1.4 Å as default.

(eq3) ∆Gbind = (Eint + Eele + Evdw ) + ( ∆GPB/GB + ∆GSA ) – T∆S

M AN U

The nonpolar solvation contribution is determined on the basis of a solvent-accessible surface area SASA solving equation 4 (eq4) where γ is a surface tension parameter and b a parameterized value. Default values were used for γ and b terms.

(eq4) ∆GSA = γ SASA + b

TE D

The final binding energy is often determined by considering the conformational entropy change expressed as (– T∆S) where T is the absolute temperature and S the entropy of the molecule. Entropy of the molecule accounts for the loss of translational, rotational and conformational degrees of free-

EP

dom of ligand upon binding and it is critical for determining the absolute binding free energies when the ligand binding is driven by large favourable entropy changes [3,17,82]. However, it has

AC C

already been reported that inclusion of entropy in calculations highly increases the computational cost, but did not always improve the accuracy [83]. It has also been shown that better correlations can be achieved upon omission of the entropy term from binding free energy calculations in case of congeneric series [25,83]. As we will be also comparing the relative binding free energies of a series of similar compounds, we did not include the entropic contribution in our calculations and only enthalpy values (∆H) were computed to see whether a correlation could be achieved between binding enthalpies and biological data (IC50 values were converted into pIC50 values; pIC50 = -logIC50) (eq5).

ACCEPTED MANUSCRIPT (eq5) ∆H = (EMM ) + ( ∆Gsol ) In addition, the estimation of conformational entropy change via the number of rotatable bonds in the ligand has been tested in docking and binding free energy calculations [22,23,84-88]. Promising

RI PT

results have been obtained implying that this measure may enhance the correlation with biological data [22,86]. The conformational entropy component was approximated from the number of rotatable bonds (NRot) of ligands (eq6) and combined with the enthalpy derived from MM-GBSA method and this sum will be referred to as MM-GBSA values.

1.6. Principal Component Analysis

M AN U

SC

(eq6) T∆SNRot = number of rotatable bonds x 1.0 kcal/mol

An analysis of the chemical space of the studied inhibitors was done using Principal Component Analysis (PCA) [89] implemented in MOE [57]. Molecular descriptors for the molecules were generated using the program MOE. A variety of 1D descriptors including the number of rotatable

TE D

single bonds (brotN), number of rings (rings), number of heavy atoms (a_heavy), number of Hbond acceptors (a_acc), number of H-bond donors (a_don), molecular weight (MW), topological polar surface area (TPSA), and log octanol/water partition coefficient (logP(o/w)), Lipinski

EP

violation count (lip_violation) were used. In order to select the important molecular descriptors that

AC C

will be used in PCA Method, a bivariate contingency analysis was performed for each descriptor and the activity data using the QuaSAR-Contingency module implemented in MOE. The selected descriptors (MW, TPSA, logP(o/w), rings, brotN, a_heavy and a_acc) obtained from the QuaSARContingency calculatuion, was further reduced by linear transformation of the data by PCA. The first three principal components obtained from PCA (PCA1, PCA2 and PCA3) were then used to plot the chemical space of the thieno[3,2-d]pyrimidine-6-carboxamide derivatives. The X axis is set to PCA1, Y to PCA2 and Z to PCA3. Each point in the three dimensional (3D) scatter plot corresponds to a molecule (Figure S6 Supporting Information).

ACCEPTED MANUSCRIPT

2. RESULTS and DISCUSSION 2.1. Analysis of Sirt3-Inhibitor Complexes

RI PT

The superimposition of the three X-ray structures of Sirt3 in complex with the reported inhibitors revealed that the inhibitors bind in a similar fashion in the catalytic binding site of Sirt3 [43]. Due to their high structural similarity this is something to be expected. The availableSirt3-11c, Sirt3-28 and Sirt3-31 complex structures unveil that the three crystallized inhibitors bind in the same way in

SC

the catalytic core pocket, occupying both acetyl lysine substrate channel as well as the nicotinamide C-pocket (which overlaps with NAD+ cofactor binding channel). A closer look to the binding

M AN U

interactions reveals that, the carboxamide group of the thieno[3,2-d]pyrimidine-6-carboxamide core mimics the corresponding group of nicotinamide (NCA) [61]. The amide group makes essential hydrogen bond interactions with the conserved residues Ile230 and Asp231 while a water-mediated H-bond interaction is formed with Ala146 and Ile154 in the C pocket. The remaining part of the

TE D

inhibitor structure adopts an extended conformation in the acetyl lysine channel benefitting from additional hydrophilic and hydrophobic contacts depending on the substituents it bears. To test whether GLIDE is able to correctly reproduce the correct binding mode of the co-crystallized

EP

inhibitors, we first re-docked the compounds 11c, 28 and 31 in the corresponding crystal structure. One conserved structural bridging water molecule (W531, W502 and W522 present in Sirt3-11c,

AC C

Sirt3-28 and Sirt3-31 complex structures, respectively) was maintained in the binding site for the docking study. The applied docking protocol used in GLIDE was able to correctly reproduce the binding mode of the inhibitor in the corresponding structure with low RMSD values (Table 1, Figure 1). To test which Sirt3 crystal structure represents the best conformation for the docking study, we carried out a cross-docking (Table 1). Sirt3-11c and Sirt3-31 structures performed well yielding an RMSD below 1.1 Å for all docked inhibitors. Using the Sirt3-28 structure GLIDE failed to correctly predict the interactions of inhibitor 11c in the substrate channel (data not shown) and gave an RMSD of 4.35 Å.

ACCEPTED MANUSCRIPT In the next step, we docked all 33 inhibitors into both Sirt3-11c and Sirt3-31 protein structures to deduce whether the docking approach is able to reproduce the correct binding mode for all compounds. All compounds could be docked into the protein structure Sirt3-31 giving the correct binding mode for the core part of the inhibitors as seen in the Sirt3 crystal structures. Meanwhile,

RI PT

compounds 17 and 39 were not correctly docked in the crystal structure Sirt3-11c (data not shown). In summary, Sirt3-31 crystal structure gave overall the best docking results and the corresponding protein-inhibitor complexes were therefore chosen for the BFE calculations (for details see

SC

Materials and Methods part).

In addition, the most potent Sirt3 inhibitor, 17, was selected to examine the stability of the

M AN U

generated complexes by means of a 20 ns MD simulations using AMBER 12. The RMSD plot showed that the complex was stable during the MD run and the RMSD values calculated for both the backbone heavy atoms of Sirt3 pocket residues and the heavy atoms of the inhibitor showed small fluctuations between 1.5 – 2.5 Å during the 20 ns of the MD simulation demonstrating the

TE D

stability of the Sirt3-inhibitor complex (see Supporting Information Figure S1).

2.2. Analysis of Sirt1-inhibitor complexes

EP

Superimposition of available Sirt1 crystal structures with Sirt3-31 showed that the Sirt1 structure complexed with NAD+ and the inhibitor Ex527 which represents an inhibited conformation (PDB

AC C

ID; 4I5I) is very similar (RMSD value of backbone atoms 1.03 Å) (Figure 2, and Supporting Information Figure S2).

In addition, the further visual inspection of the Sirt1-NAD+-Ex527 crystal structure revealed that the inhibitor Ex527 is showing a similar binding mode as the reported Sirt3 inhibitors 11c, 28 and 31 (hydrogen bond interactions with the conserved residues Ile347 and Asp348 and the bridging water in the catalytic binding site). These findings prompted us to choose the Sirt1-NAD+-Ex527 crystal structure for the docking study. The conformation of the residues lining the binding pocket observed in Sirt3-31 and Sirt1-NAD+-Ex527 is similar except for the conserved residue Phe157

ACCEPTED MANUSCRIPT (Sirt3 numbering). In the Sirt3-31 crystal structure, Phe157 is making π –π interactions with the thienopyrimidine core of the inhibitor and thus forcing the compound to get a bended conformation and forming hydrogen bond interactions with the conserved residues of the binding pocket. To maintain this critical role of the conserved Phe, the side-chain conformation of Phe273 in Sirt1-

RI PT

NAD+-Ex527 was changed using the rotamer library implemented in MOE. As consequence, the inhibitor 31 could be docked into Sirt1 resulting in the same binding mode as observed for Sirt3. Furthermore, we applied a 1 ns MD simulation of the Sirt1-31 complex structure to relax the

SC

protein conformation. The Sirt1 structure which was obtained after the 1ns MDs was then used for the further docking and binding free energy calculations. Only 31 inhibitors out of 33 could be

M AN U

correctly docked into the Sirt1 protein. No poses could be obtained for compound 17 and 42, respectively. The most potent Sirt1 inhibitor, 11c, was selected to examine the stability of the generated complexes by means of a 20 ns MD simulations using AMBER 12. The RMSD plot showed that the complex was stable during the MD simulation and the RMSD values calculated for both the backbone heavy atoms of the pocket residues and the heavy atoms of the inhibitor

TE D

remained stable around 1.5 Å (see Supporting Information Figure S3).

EP

2.3. Analysis of Sirt2-Inhibitor Complexes

The same approach as used for Sirt1 was applied to find the most appropriate Sirt2 conformation

AC C

out of the already published Sirt2 crystal structures. The co-substrate analog, ADPR, bound Sirt2 crystal structure was found to be most similar to Sirt3-31 (Figure 3). An RMSD value of 1.06 Ǻ (backbone atoms) was calculated (see Supporting Information Figure S4). Superimposition of Sirt2 and Sirt3 structures revealed that the conserved bridging water is also present in the Sirt2-ADPR structure (A chain; W2024) and shows an identical position as the water molecule (W522) in the Sirt3-31 complex. Thus, we kept the structural water, W2024, at the same position in the Sirt2-ADPR structure for further modelling studies. A suitable conformation of the conserved residue Phe96 was obtained using the rotamer library implemented in MOE as

ACCEPTED MANUSCRIPT mentioned before. Subsequently, the Sirt2 protein structure in complex with the potent and large inhibitor 15a was subjected to a short MD simulation. The docking studies were carried out using the MD-derived Sirt2 protein conformation and all 33 compounds could be correctly docked. Protein-ligand complex generation and subsequent BFE calculation were performed for each ligand

RI PT

as described before (see the Materials and Methods part). In addition the stability of the Sirt2inhibitor complex was analyzed by running a long MD simulation (20 ns). The most potent Sirt2 inhibitor among the series (31) was used for this purpose. The RMSD plot of the protein and

SC

inhibitor 31 showed that the complex was stable during the complete MD run and the RMSD values calculated for both backbone heavy atoms of Sirt2 pocket residues and the heavy atoms of the

M AN U

inhibitor showed small fluctuations between 1.5 – 2.5 Å during the 20 ns production MD period demonstrating the stability of the Sirt2-inhibitor binding (see Supporting Information Figure S5).

2.4. Docking results

TE D

Docking results of the 33 studied inhibitors showed that the compounds display a similar binding mode in all Sirtuin binding pockets as seen in the inhibitor bound crystal structures of Sirt3. Figures 4 and 5 show the interaction of the most potent inhibitor of this series at Sirt1 and Sirt2.

EP

The compounds are forming H-bonds with the conserved residues Ile347/169 and Asp348/170 (Sirt1 and Sirt2 numbering, respectively) and are also benefiting from additional interactions

AC C

through the structural water present in the binding site. While the binding modes of the studied inhibitors could be correctly predicted, no correlation could be attained between the docking scores (GLDE-XP) and the pIC50 values (R2 = 0.26 for Sirt1, R2 = 0.18 for Sirt2 and R2 = 0.16 for Sirt3, Figure 6). The docking results are summarized in Table 2, 3 and 4. Often, docking scores show a poor correlation with the experimental inhibition or binding data. On the other hand, we also realised that the obtained GLIDE-XP scores give a hint for classifying the compounds into active and inactive inhibitors (compounds showing IC50 values >50 µM). Box plots are generated showing the distribution of GLIDE-XP scores obtained for active and

ACCEPTED MANUSCRIPT inactive compounds for Sirt1, Sirt2 and Sirt3 (Figure 7). It can be seen in the box plots that GLIDE-XP scores obtained for inactive compounds are considerably lower than the active compounds in general for all sirtuin subtypes under study. A median value of -7.85/-6.57, -9.58/-

Sirt3, respectively.

SC

2.5. Rescoring of docking poses using MM-PB(GB)SA

RI PT

7.54 and -11.30/-8.62 for the active/inactive compounds was acquired in case of Sirt1, Sirt2 and

Since docking scores often show no correlation to biological in vitro data, more sophisticated

M AN U

binding free energy calculation methods are sometimes applied to yield a better performance [18,23,83,90]. Rescoring by applying MM-GBSA or MM-PBSA can be performed either using single protein-inhibitor complexes or by considering several frames derived from MD simulations [23,91]. In the current study we tested several setups for binding free energy calculation applying different MM-PB(GB)SA solvation models (MM-PBSA and MM-GBSA method will be referred

TE D

according to its ipb or igb value PB-1 (ipb=1), GB-1 (igb=1), GB-2 (igb=2), GB-5 (igb=5), GB-8 (igb=8)). Protein-ligand complexes were prepared for all compounds for which a reported IC50 value was available. Thus, in total 27 complexes for Sirt1, 30 complexes for Sirt2 and 28 complexes

EP

for Sirt3 were submitted to BFE calculations. The poses were first subjected to an energy minimization step using an explicit solvent model and then a short equilibration MD run of 100 ps.

AC C

The motivation behind this is that by introducing an MD step the conformational space obtained from one single docking pose will be enlarged and thus, the probability of getting complexes of ligands showing a better fit to the protein can be increased. Moreover, averaging the binding energies retrieved from a small collection of protein-inhibitor complexes can help to alleviate the error that comes from the energy terms which are sensitive to the exact position of the ligand in the binding (e.g. the vdW term). In addition, we tested three different setups to get the final binding free energy values. In the first setup, we used the protein-inhibitor conformation after energy minimization in explicit water to calculate the enthalpy values. In the second setup, we used the

ACCEPTED MANUSCRIPT averaged enthalpy value calculated for 10 conformations obtained from the 100 ps MD run and in the third setup the conformation obtained after 100 ps MD was used. The performance of the BFE methods was evaluated by computing the correlation between

RI PT

biological activity data and enthalpy values derived from each method. The data are summarized in Table S1, S2 and S3 (see Supporting Information) and the correlation plots are presented in Figure 8. Besides the correlation coefficient R2 the models were also tested using LOO (Leave One Out) cross-validation. Q2LOO and the XZ-scores were calculated in addition. In this application, the XZ-

SC

score represents the absolute difference between the predicted value of the model under a leaveone-out cross validation scheme and the activity data, divided by the square root of the mean square

M AN U

error of the data set. Molecules with a calculated XZ-score of 2.0 or higher were defined as outliers and removed from the dataset to further refine the model.

The statistical analyses obtained for Sirt3 are summarized in Table 5. The results for the Sirt3 inhibitors indicate that after removing the outliers (11c, 23, 31 and 38), the GB-1 model using 10

TE D

MD frames (GB1-md) gives the best correlation overall (R2 = 0.60, RMSE = 0.50, q2LOO = 0.54, n = 24). The second best performance was obtained by using GB-8 model either by using the final conformation acquired after minimization (GB8-min) (R2 = 0.57, RMSE = 0.51, q2LOO = 0.47, n =

EP

24) or using 10 MD frames (GB8-md) (R2 = 0.57, RMSE = 0.51, q2LOO = 0.48, n = 24). GB-1 and GB-8 models were found to get better results compared to the other solvation models. Therefore,

AC C

GB-1 and GB-8 models were chosen for further binding free energy calculations using the Sirt1 and Sirt2 inhibitors as well. Binding free energies and statistics obtained from GB-1 and GB-8 calculations of Sirt1 and Sirt2 complexes are shown in Table S4 and Table S5, (see Supporting Information) and Table 6 and Table 7, respectively. Binding free energies obtained from the GB1-md model for Sirt1, Sirt2 and Sirt3 are presented in Table 2, 3 and 4, respectively. The results show that GB1-md model also performed best for Sirt1 (R2 = 0.79, RMSE = 0.41, q2LOO = 0.76, n = 24) and Sirt2 inhibitors (R2 = 0.76, RMSE = 0.51, q2LOO = 0.72, n = 27). A significant improvement can be observed compared to the scoring-based

ACCEPTED MANUSCRIPT models (GLIDE-XP: R2 = 0.26 for Sirt1, R2 = 0.18 for Sirt2 and R2 = 0.16 for Sirt3). A rough discrimination between active and inactive compounds can be obtained by using the calculated enthalpy values. As shown in Figure 9, enthalpy scores retrieved for active compounds are

RI PT

more favourable than those for inactive compounds. A median value of 53.79/41.22, 60.00/43.60 and 56.36/43.98 kcal/mol for the active/inactive compounds was acquired in case of Sirt1, Sirt2 and Sirt3, respectively. We observed that the inactive compounds 14 and 22 (IC50 values determined for compounds 14 and 22 against Sirt1 and Sirt3 are both >50µM) could be easily discriminated from

SC

actives according to their enthalpy scores in case of Sirt1 and Sirt3 and ranked at the bottom in case of Sirt2 (IC50 values determined for compounds 14 and 22 against Sirt2 are >50 µM and 49µM, re-

M AN U

spectively).

The weakly active inhibitors of Sirt1 (23), Sirt2 (41), and Sirt3 (23 and 38) were pointed out as outliers in the cross-validation models. One reason can be the low activity (several orders of magnitude less active) of these inhibitors which is not reflected by the derived BFE models. A

TE D

reasonable explanation why the potent inhibitors of Sirt1 (19 and 31), Sirt2 (30 and 31), and Sirt3 (11c and 31) were indicated as outliers still has to be elucidated. As we already mentioned before, the energy terms are sensitive to the exact position of the ligand in

EP

the binding pocket (e.g. the vdW term) and even a slight difference in conformations of the ligand may result in lower or higher enthalpy scores. Moreover, losing essential binding interactions

AC C

and/or forming misleading interactions in the binding pocket can deteriorate the final enthalpy score entirely. The effectiveness of binding free energy values to predict the correct protein-ligand binding has already been assessed in literature and examples of promising results are reported [87,92,93]. It was shown that the protein-ligand complexes with the lowest binding energy had the closest match to the X-ray structure rather than the top-ranked docking pose [21]. We sought to determine whether we could get better MM-GBSA values and thus, better binding poses by postprocessing all the docked poses produced by GLIDE docking for as the identified outliers. In case of Sirt1, the top-ranked poses used for outliers 19 and 31 in MM-GBSA calculations were

ACCEPTED MANUSCRIPT showing all the H-bond interactions between thieno[3,2-d]pyrimidine-6-carboxamide core and the conserved residues Ile347 and Asp348 and the water molecule needed for the inhibitory activity, however, the H-bond between backbone C=O of Val292 and the NH of group R3 (see Chart1) seen in the crystal structures was not present [43]. We observed the same findings upon analysing the

RI PT

protein-ligand complexes with the lowest ∆H values for outliers 19 and 31. Yet, these complexes were sharing a more similar conformation of the piperidine moiety with the crystal structure and this might be the reason of the increase in the enthalpy score (Table 8).

SC

The same approach as used for Sirt1 was applied to Sirt2 outliers (30 and 31). GLIDE-XP returned only one docking pose for inhibitor 31 and therefore, no comparison could be made. For compound

M AN U

30, a better enthalpy score could be obtained upon including all the docking poses in MM-GBSA calculations (Table8). While the substantial interactions needed for the sirtuin inhibitory activity discussed above were correctly assigned for 30, the H-bond interaction with Val233 could not be regained with the top-ranked docking pose. On the other hand, the Sirt2-30 complex with the lowest

TE D

∆H value was able to get the H-bond interaction with the residue Val233 (data not shown) and hence improved enthalpy values (Table8).

In general, there is considerable increase in the absolute enthalpy values and improvement in the

EP

binding poses upon considering all docking poses and ranking them by ∆H values rather than focusing exclusively on the top-ranked pose in case of Sirt1 and Sirt2. The results suggest that the

AC C

loss of H-bond interaction with the Val residue and/or slight differences in conformation of piperidine moiety was the reason behind the unfavourable enthalpy values retrieved for the outliers. Nonetheless, the results are controversial in case of Sirt3. We realised that there is also an improvement in enthalpy values upon taking into consideration all docking poses generated for the outliers of Sirt3, compound 11c and 31 (Table 8). To explore the competence of each method, we superimposed the top-ranked poses acquired with GLIDE scoring and MM-GBSA method with the crystals of Sirt3-11c and Sirt3-31 (Figures 10 and 11) in MOE. As depicted in Figure 10, compound 11c, could be correctly docked in the binding pocket getting

ACCEPTED MANUSCRIPT all the fundemental interactions coming from the thieno[3,2-d]pyrimidine-6-carboxamide core. On the other hand, a closer look to the top-ranked docking pose revealed that while the H-bond between the distal ethylamide substituent on the 2-thiophene group of 11c and Glu296 could be correctly predicted, the H-bond between the aryl amide NH and Val292 in the substrate channel

RI PT

could not be regained. The absence of H-bond interaction between the NH donor and the Val292 might explain the poor ∆H score acquired for compound 11c and thus, being indicated as an outlier. On the contrary, the ethylamide substituted 2-thiophene moiety did not show the extended

SC

conformation seen in the crystal structure in case of the lowest binding energy pose and formed a H-bond with the residue His248 instead (Figure 10). In conclusion, replacing the loss of H-bond

M AN U

interactions expected with the Val292 and Glu296 residues with a wrongly predicted H-bond with His248 residue was well tolerated, resulting in modest improvement in the enthalpy score (Table 8). Taking the case of 11c, the hypothesis “the lowest binding energy pose likely to have the closest match to the X-ray structure rather than the top-ranked docking pose” appears to be insufficient.

TE D

Superimposition of the co-crystallized ligand 31 with the top-ranked docking pose revealed that the ethyl linker situated between the piperidine moiety and the sulphonamide group is twisted and a different conformation was attained by the sulphonamide group in the top-ranked docking pose with

EP

an RMSD of 1.00 Å (Figure 11). On the other hand, a more similar conformation as in the crystal structure could be obtained for the ethyl linker in case of the lowest binding energy pose getting an

AC C

RMSD of 0.89 Å (Figure 11). Since all the other significant binding interactions beside the ethyl linker and sulphonamide group could be correctly predicted in both cases, we suppose that the improvement in enthalpy score comes from the better positioning of ethyl linker in the lowest binding An increase in the enthalpy values have been shown for most of the cases when the lowest binding energy pose was used instead of the top-ranked docking pose. To test whether there is a better correlation between experimental pIC50 values and ∆H values derived from lowest binding free energy poses, we carried out MM-GBSA calculations using all poses obtained from docking into the Sirt3 structure. Subsequently, we took the pose with the lowest binding free energy for the statistical

ACCEPTED MANUSCRIPT analysis. Stastical results are summarized in Table S6 (see Supporting Information). We found a slight increase of R2 from 0.60 to 0.62 (n=24) and q2LOO from 0.54 to 0.56 (n=24) upon using the lowest binding energy pose for correlation instead of the top-ranked docking pose. Moreover, as already discussed the lowest binding energy pose is not always close to the crystal

RI PT

structure compared to the top-ranked docking pose. In addition, we showed that with less computing time a good correlation can be already achieved between the biological data and the binding

SC

free energies using the top-ranked docking pose.

2.6. Applicability Domain

M AN U

In order to define the chemical space of the inhibitor data set, a Principal Components Analysis (PCA) [89] was performed. PCA represents a method for describing the Applicability Domain (AD) for a given Quantitative Structure-Activity Relationships (QSAR) model [94,95]. The AD determines the limitations of the model with respect to the chemical space of the studied compounds [94]. In the present study, the AD is defined by three principal components (PCA1, PCA2 and

TE D

PCA3) calculated from simple descriptors (MW, TPSA, logP(o/w), rings, brotN, a_heavy and a_acc) of the studied inhibitors having an exact IC50 value for Sirt3. To test the AD we took the Sirt3 results as data set [96]. The PCA plot of the inhibitor data set for Sirt3 is shown in Figure S6

EP

(see Supporting Information).

Five different groups of molecules were observed in the PCA plot, named as G1, G2, G3, G4 and

AC C

G5. Group G1 consists of most of the molecules in the inhibitor data set (15a, 16a, 20, 23, 24, 25, 26, 27, 28, 29, 30, 31, 32, 36, 37, 38, 39). These compounds were shown to inhibit Sirt3 with an IC50 in the range of 0.0072 to 4.4 µM. Group G2 consists of only compound 21, which is the least active compound with an exact IC50 value of 8.6 µM. It does not contain R3 group (see Chart1) and therefore, it has a smaller size compared to the rest of the compounds present in G1, G3, G4 and G5. Group G3 consists of three compounds; 33, 34 and 35. These compounds were shown to inhibit Sirt3 with an IC50 in the range of 0.028 to 0.061 µM. Group G4 consists of only compound

ACCEPTED MANUSCRIPT 19, which inhibits Sirt3 with an IC50 value of 0.0076 µM. Group G5 consists of compounds 11a, 11b, 11c, 11d, 17, and 18. These compounds were shown to inhibit Sirt3 with an IC50 in the range of 0.0032 to 0.024 µM. In MM-GBSA method compound 11c, 23, 38 and , 31 were predicted as outliers which reside in

RI PT

G1 and G5 in the PCA plot, respectively. Since the outliers are located within the observed five groups in the PCA plot, the reason why these compounds were found to be as outliers in the MMGBSA model is not related to their physicochemical properties.

SC

Even though, the inhibitor data set is made of congeneric series of thieno[3,2-d]pyrimidine-6carboxamide derivatives, we show by calculating the PCA plot that compounds cover a wide range

M AN U

of different physicochemical properties which also implies that the AD is not too restricted. These findings strengthen the hypothesis that small structural changes observed in the inhibitor structure can on the one site significantly affect the biological activity. In addition, the docking and binding free energy calculations of structurally similar compounds which are equally active (e.g. 11a-11d on

lation.

TE D

Sirt3, might result in a misprediction due to the loss of important interactions during the MD simu-

EP

2.7. Evaluation of entropy contribution

We also investigated the influence of including the entropy term estimated by using the number of

AC C

rotatable bonds of the compounds (see Table 9 for a summary of statistics). We observed a slight decrease in the absolute values of the correlation coefficient and cross-correlation coefficient values for all Sirtuin subtypes upon including the entropy term. In case of Sirt1, R2 decreased from 0.79 to 0.75 and q2LOO decreased from 0.76 to 0.71 (n=24). In case of Sirt2 and Sirt3, there is also a considerable decrease of R2 from 0.76 to 0.72 and for q2LOO from 0.72 to 0.66 (Sirt2, n=28) and R2 from 0.60 to 0.54 and for q2LOO from 0.54 to 0.47 (Sirt3, n=25). However, this might be due to the increase in the final dataset from 27 compounds to 28 for Sirt2 (compound 30 which was found to be an outlier in the enthalpy model is included in this model) and from 24 compounds to 25 for

ACCEPTED MANUSCRIPT Sirt3 (compound 11c which was found to be an outlier in the enthalpy model is included in this model) that were used for the correlation.

RI PT

3. CONCLUSION In this study, we explored the performance of different implicit solvation models, as post-processing tool for rescoring docking poses derived for a dataset of thieno[3,2-d]pyrimidine-6-carboxamide

SC

derivatives tested as Sirt1-3 inhibitors. Various setups were combined to produce the protein-ligand complexes used for the estimation of the binding free energies. Since there were no Sirt1 or Sirt2

M AN U

crystal structures available for the inhibitors under study, in the first step, we generated proteininhibitor modes for each enzyme. We showed that the correct bound conformation of the Sirt isoform could be employed in a very fast and accurate way by applying a short MD simulation. In the second part of the study, we showed that a moderate discrimination between inactive and active compounds could be achieved by using the GLIDE-XP docking scores. However, only a poor

TE D

or no correlation between the docking scores and the experimental data of the active compounds was observed. To establish a better performing model, we took the docking poses of the inhibitors in the relevant sirtuin structures and applied energy minimization and a subsequent short MD

EP

simulation. Two models either based on MM-GBSA enthalpy values or entropy included MMGBSA values were generated to test whether the derived binding free energies could correctly rank

AC C

the compounds according to their determined biological activity. Testing different MM-GBSA models resulted in a final setup that showed high correlation and leave-one-out cross correlation coefficients. The best overall performance was observed using the GB-1 solvation model and ten frames derived from a short MD simulation. We showed that promising results were obtained using the fairly new MM-GBSA GB-8 solvation model [62] which is not well studied in literature until now. Finally, the effect of including the entropy contribution to the MM-GBSA method was assessed. Instead of calculating the entropy term with computationally demanding normal mode analysis, we took a simple approach where just the number of rotatable bonds is taken into account.

ACCEPTED MANUSCRIPT Addition of the entropy term did not improve the overall performance of the MM-GBSA. Whether more demanding normal mode analysis on the derived Sirt-inhibitor will result in better performing models has to be tested in future work.

RI PT

The derived MM-GBSA model represents a time-efficient post-docking filter also for larger data sets. The CPU time required for the calculation of one inhibitor including energy minimization of the solvated protein-ligand complex, equilibration MD run and binding free energy calculation is around ~ 10 minutes. Therefore, the docking poses of thousands of compounds can be rescored in a

SC

reasonable amount of time. Second, with respect to the MM-GBSA model, an overall better performance was found compared with the GLIDE-XP docking model in terms of discrimination of

M AN U

active compounds from inactive compounds. Most important, rescoring of docked poses with MMGBSA method showed a better agreement with the biological data than using docking scores.

ACKNOWLEDGEMENT

TE D

Authors would like to give acknowledgements to the Institute of Pharmaceutical Chemistry, Martin-

AC C

EP

Luther University Halle-Wittenberg, Halle (Saale), Germany.

ACCEPTED MANUSCRIPT REFERENCES

AC C

EP

TE D

M AN U

SC

RI PT

[1] Sherman, W.; Day, T.; Jacobson, M. P.; Friesner, R. A.; Farid, R., Novel procedure for modeling ligand/receptor induced fit effects. J Med Chem 2006, 49 (2), 534-53. [2] Warren, G. L.; Andrews, C. W.; Capelli, A. M.; Clarke, B.; LaLonde, J.; Lambert, M. H.; Lindvall, M.; Nevins, N.; Semus, S. F.; Senger, S.; Tedesco, G.; Wall, I. D.; Woolven, J. M.; Peishoff, C. E.; Head, M. S., A critical assessment of docking programs and scoring functions. J Med Chem 2006, 49 (20), 5912-31. [3] Gilson, M. K.; Zhou, H. X., Calculation of protein-ligand binding affinities. Annu Rev Biophys Biomol Struct 2007, 36, 21-42. [4] Moitessier, N.; Englebienne, P.; Lee, D.; Lawandi, J.; Corbeil, C. R., Towards the development of universal, fast and highly accurate docking/scoring methods: a long way to go. Br J Pharmacol 2008, 153 Suppl 1, S7-26. [5] Vaque, M. A., Anna; Blade, Cinta; Salvado, M. J.; Blay, Mayte; Fernandez-Larrea, Juan; Arola, Lluis, Pujadas, Gerard, Protein-ligand Docking: A Review of Recent Advances and Future Perspectives. Current Pharmaceutical Analysis 2008, 4, 1-19. [6] Armen, R. S.; Chen, J.; Brooks, C. L., 3rd, An Evaluation of Explicit Receptor Flexibility in Molecular Docking Using Molecular Dynamics and Torsion Angle Molecular Dynamics. J Chem Theory Comput 2009, 5 (10), 2909-2923. [7] Soulere, L., Toward docking-based virtual screening for discovering antitubulin agents by targeting taxane and colchicine binding sites. ChemMedChem 2009, 4 (2), 161-3. [8] Durrant, J. D.; McCammon, J. A., Molecular dynamics simulations and drug discovery. BMC Biol 2011, 9, 71. [9] Wichapong, K.; Lindner, M.; Pianwanit, S.; Kokpol, S.; Sippl, W., Receptor-based 3DQSAR studies of checkpoint Wee1 kinase inhibitors. Eur J Med Chem 2009, 44 (4), 1383-95. [10] Yuriev, E.; Agostino, M.; Ramsland, P. A., Challenges and advances in computational docking: 2009 in review. J Mol Recognit 2011, 24 (2), 149-64. [11] Aqvist, J.; Medina, C.; Samuelsson, J. E., A new method for predicting binding affinity in computer-aided drug design. Protein Eng 1994, 7 (3), 385-91. [12] Srinivasan, J.; Cheatham, T. E.; Cieplak, P.; Kollman, P. A.; Case, D. A., Continuum Solvent Studies of the Stability of DNA, RNA, and Phosphoramidate−DNA Helices. J Am Chem Soc 1998, 120 (37), 9401-9409. [13] Kirkwood, J. G., Statistical Mechanics of Fluid Mixtures. The Journal of Chemical Physics 1935, 3 (5), 300-313. [14] Straatsma, T. P.; McCammon, J. A., Multiconfiguration thermodynamic integration. The Journal of Chemical Physics 1991, 95 (2), 1175-1188. [15] Archontis, G.; Watson, K. A.; Xie, Q.; Andreou, G.; Chrysina, E. D.; Zographos, S. E.; Oikonomakos, N. G.; Karplus, M., Glycogen phosphorylase inhibitors: A free energy perturbation analysis of glucopyranose spirohydantoin analogues. Proteins: Structure, Function, and Bioinformatics 2005, 61 (4), 984-998. [16] Pearlman, D. A., Evaluating the Molecular Mechanics Poisson−Boltzmann Surface Area Free Energy Method Using a Congeneric Series of Ligands to p38 MAP Kinase. J Med Chem 2005, 48 (24), 7796-7807. [17] Homeyer, N.; Gohlke, H., Free Energy Calculations by the Molecular Mechanics Poisson−Boltzmann Surface Area Method. Molecular Informatics 2012, 31 (2), 114-122. [18] Kuhn, B.; Gerber, P.; Schulz-Gasch, T.; Stahl, M., Validation and use of the MM-PBSA approach for drug discovery. J Med Chem 2005, 48 (12), 4040-8. [19] Rastelli, G.; Del Rio, A.; Degliesposti, G.; Sgobba, M., Fast and accurate predictions of binding free energies using MM-PBSA and MM-GBSA. J Comput Chem 2010, 31 (4), 797-810. [20] Wichapong, K.; Lawson, M.; Pianwanit, S.; Kokpol, S.; Sippl, W., Postprocessing of protein-ligand docking poses using linear response MM-PB/SA: application to Wee1 kinase inhibitors. J Chem Inf Model 2010, 50 (9), 1574-88.

ACCEPTED MANUSCRIPT

AC C

EP

TE D

M AN U

SC

RI PT

[21] Lindstrom, A.; Edvinsson, L.; Johansson, A.; Andersson, C. D.; Andersson, I. E.; Raubacher, F.; Linusson, A., Postprocessing of docked protein-ligand complexes using implicit solvation models. J Chem Inf Model 2011, 51 (2), 267-82. [22] Slynko, I.; Scharfe, M.; Rumpf, T.; Eib, J.; Metzger, E.; Schule, R.; Jung, M.; Sippl, W., Virtual screening of PRK1 inhibitors: ensemble docking, rescoring using binding free energy calculation and QSAR model development. J Chem Inf Model 2014, 54 (1), 138-50. [23] Wichapong, K.; Rohe, A.; Platzer, C.; Slynko, I.; Erdmann, F.; Schmidt, M.; Sippl, W., Application of docking and QM/MM-GBSA rescoring to screen for novel Myt1 kinase inhibitors. J Chem Inf Model 2014, 54 (3), 881-93. [24] Rastelli, G.; Degliesposti, G.; Del Rio, A.; Sgobba, M., Binding Estimation after Refinement, a New Automated Procedure for the Refinement and Rescoring of Docked Ligands in Virtual Screening. Chemical Biology & Drug Design 2009, 73 (3), 283-286. [25] Uciechowska, U.; Schemies, J.; Scharfe, M.; Lawson, M.; Wichapong, K.; Jung, M.; Sippl, W., Binding free energy calculations and biological testing of novel thiobarbiturates as inhibitors of the human NAD(+) dependent histone deacetylase Sirt2. Medchemcomm 2012, 3 (2), 167-173. [26] Jiang, H.; Khan, S.; Wang, Y.; Charron, G.; He, B.; Sebastian, C.; Du, J.; Kim, R.; Ge, E.; Mostoslavsky, R.; Hang, H. C.; Hao, Q.; Lin, H., SIRT6 regulates TNF-alpha secretion through hydrolysis of long-chain fatty acyl lysine. Nature 2013, 496 (7443), 110-3. [27] Du, J.; Zhou, Y.; Su, X.; Yu, J. J.; Khan, S.; Jiang, H.; Kim, J.; Woo, J.; Kim, J. H.; Choi, B. H.; He, B.; Chen, W.; Zhang, S.; Cerione, R. A.; Auwerx, J.; Hao, Q.; Lin, H., Sirt5 Is a NADDependent Protein Lysine Demalonylase and Desuccinylase. Science 2011, 334 (6057), 806-809. [28] Ahuja, N.; Schwer, B.; Carobbio, S.; Waltregny, D.; North, B. J.; Castronovo, V.; Maechler, P.; Verdin, E., Regulation of insulin secretion by SIRT4, a mitochondrial ADP-ribosyltransferase. J Biol Chem 2007, 282 (46), 33583-92. [29] North, B. J.; Verdin, E., Sirtuins: Sir2-related NAD-dependent protein deacetylases. Genome Biol 2004, 5 (5), 224. [30] Flick, F.; Lüscher, B., Regulation of sirtuin function by posttranslational modifications. Frontiers in Pharmacology 2012, 3. [31] Hiratsuka, M.; Inoue, T.; Toda, T.; Kimura, N.; Shirayoshi, Y.; Kamitani, H.; Watanabe, T.; Ohama, E.; Tahimic, C. G.; Kurimasa, A.; Oshimura, M., Proteomics-based identification of differentially expressed genes in human gliomas: down-regulation of SIRT2 gene. Biochem Biophys Res Commun 2003, 309 (3), 558-66. [32] Outeiro, T. F.; Kontopoulos, E.; Altmann, S. M.; Kufareva, I.; Strathearn, K. E.; Amore, A. M.; Volk, C. B.; Maxwell, M. M.; Rochet, J. C.; McLean, P. J.; Young, A. B.; Abagyan, R.; Feany, M. B.; Hyman, B. T.; Kazantsev, A. G., Sirtuin 2 inhibitors rescue alpha-synuclein-mediated toxicity in models of Parkinson's disease. Science 2007, 317 (5837), 516-9. [33] Harting, K.; Knoll, B., SIRT2-mediated protein deacetylation: An emerging key regulator in brain physiology and pathology. Eur J Cell Biol 2010, 89 (2-3), 262-9. [34] Alhazzazi, T. Y.; Kamarajan, P.; Joo, N.; Huang, J. Y.; Verdin, E.; D'Silva, N. J.; Kapila, Y. L., Sirtuin-3 (SIRT3), a novel potential therapeutic target for oral cancer. Cancer 2011, 117 (8), 1670-8. [35] Hirschey, M. D.; Shimazu, T.; Jing, E.; Grueter, C. A.; Collins, A. M.; Aouizerat, B.; Stancakova, A.; Goetzman, E.; Lam, M. M.; Schwer, B.; Stevens, R. D.; Muehlbauer, M. J.; Kakar, S.; Bass, N. M.; Kuusisto, J.; Laakso, M.; Alt, F. W.; Newgard, C. B.; Farese, R. V., Jr.; Kahn, C. R.; Verdin, E., SIRT3 deficiency and mitochondrial protein hyperacetylation accelerate the development of the metabolic syndrome. Mol Cell 2011, 44 (2), 177-90. [36] Stunkel, W.; Campbell, R. M., Sirtuin 1 (SIRT1): the misunderstood HDAC. J Biomol Screen 2011, 16 (10), 1153-69. [37] Zhang, F.; Wang, S.; Gan, L.; Vosler, P. S.; Gao, Y.; Zigmond, M. J.; Chen, J., Protective effects and mechanisms of sirtuins in the nervous system. Prog Neurobiol 2011, 95 (3), 373-95. [38] Pagans, S.; Pedal, A.; North, B. J.; Kaehlcke, K.; Marshall, B. L.; Dorr, A.; Hetzer-Egger, C.; Henklein, P.; Frye, R.; McBurney, M. W.; Hruby, H.; Jung, M.; Verdin, E.; Ott, M., SIRT1

ACCEPTED MANUSCRIPT

AC C

EP

TE D

M AN U

SC

RI PT

regulates HIV transcription via Tat deacetylation. PLoS Biol 2005, 3 (2), e41. [39] Bruzzone, S.; Parenti, M. D.; Grozio, A.; Ballestrero, A.; Bauer, I.; Del Rio, A.; Nencioni, A., Rejuvenating sirtuins: the rise of a new family of cancer drug targets. Curr Pharm Des 2013, 19 (4), 614-23. [40] Andreoli, F.; Barbosa, A. J.; Parenti, M. D.; Del Rio, A., Modulation of epigenetic targets for anticancer therapy: clinicopathological relevance, structural data and drug discovery perspectives. Curr Pharm Des 2013, 19 (4), 578-613. [41] Finnin, M. S.; Donigian, J. R.; Pavletich, N. P., Structure of the histone deacetylase SIRT2. Nat Struct Biol 2001, 8 (7), 621-5. [42] Szczepankiewicz, B. G.; Dai, H.; Koppetsch, K. J.; Qian, D.; Jiang, F.; Mao, C.; Perni, R. B., Synthesis of carba-NAD and the structures of its ternary complexes with SIRT3 and SIRT5. J Org Chem 2012, 77 (17), 7319-29. [43] Disch, J. S.; Evindar, G.; Chiu, C. H.; Blum, C. A.; Dai, H.; Jin, L.; Schuman, E.; Lind, K. E.; Belyanskaya, S. L.; Deng, J.; Coppo, F.; Aquilani, L.; Graybill, T. L.; Cuozzo, J. W.; Lavu, S.; Mao, C.; Vlasuk, G. P.; Perni, R. B., Discovery of thieno[3,2-d]pyrimidine-6-carboxamides as potent inhibitors of SIRT1, SIRT2, and SIRT3. J Med Chem 2013, 56 (9), 3666-79. [44] Gertz, M.; Fischer, F.; Nguyen, G. T.; Lakshminarasimhan, M.; Schutkowski, M.; Weyand, M.; Steegborn, C., Ex-527 inhibits Sirtuins by exploiting their unique NAD+-dependent deacetylation mechanism. Proc Natl Acad Sci U S A 2013, 110 (30), E2772-81. [45] Moniot, S.; Schutkowski, M.; Steegborn, C., Crystal structure analysis of human Sirt2 and its ADP-ribose complex. J Struct Biol 2013, 182 (2), 136-43. [46] Nguyen, G. T.; Gertz, M.; Steegborn, C., Crystal structures of sirt3 complexes with 4'bromo-resveratrol reveal binding sites and inhibition mechanism. Chem Biol 2013, 20 (11), 137585. [47] Nguyen, G. T.; Schaefer, S.; Gertz, M.; Weyand, M.; Steegborn, C., Structures of human sirtuin 3 complexes with ADP-ribose and with carba-NAD+ and SRT1720: binding details and inhibition mechanism. Acta Crystallogr D Biol Crystallogr 2013, 69 (Pt 8), 1423-32. [48] Zhao, X.; Allison, D.; Condon, B.; Zhang, F.; Gheyi, T.; Zhang, A.; Ashok, S.; Russell, M.; MacEwan, I.; Qian, Y.; Jamison, J. A.; Luz, J. G., The 2.5 A crystal structure of the SIRT1 catalytic domain bound to nicotinamide adenine dinucleotide (NAD+) and an indole (EX527 analogue) reveals a novel mechanism of histone deacetylase inhibition. J Med Chem 2013, 56 (3), 963-9. [49] Davenport, A. M.; Huber, F. M.; Hoelz, A., Structural and functional analysis of human SIRT1. J Mol Biol 2014, 426 (3), 526-41. [50] Yamagata, K.; Goto, Y.; Nishimasu, H.; Morimoto, J.; Ishitani, R.; Dohmae, N.; Takeda, N.; Nagai, R.; Komuro, I.; Suga, H.; Nureki, O., Structural basis for potent inhibition of SIRT2 deacetylase by a macrocyclic peptide inducing dynamic structural change. Structure 2014, 22 (2), 345-52. [51] Pan, P. W.; Feldman, J. L.; Devries, M. K.; Dong, A.; Edwards, A. M.; Denu, J. M., Structure and biochemical functions of SIRT6. J Biol Chem 2011, 286 (16), 14575-87. [52] Lawson, M.; Uciechowska, U.; Schemies, J.; Rumpf, T.; Jung, M.; Sippl, W., Inhibitors to understand molecular mechanisms of NAD(+)-dependent deacetylases (sirtuins). Biochim Biophys Acta 2010, 1799 (10-12), 726-39. [53] Napper, A. D.; Hixon, J.; McDonagh, T.; Keavey, K.; Pons, J. F.; Barker, J.; Yau, W. T.; Amouzegh, P.; Flegg, A.; Hamelin, E.; Thomas, R. J.; Kates, M.; Jones, S.; Navia, M. A.; Saunders, J. O.; DiStefano, P. S.; Curtis, R., Discovery of indoles as potent and selective inhibitors of the deacetylase SIRT1. J Med Chem 2005, 48 (25), 8045-54. [54] Suzuki, T.; Khan, M. N.; Sawada, H.; Imai, E.; Itoh, Y.; Yamatsuta, K.; Tokuda, N.; Takeuchi, J.; Seko, T.; Nakagawa, H.; Miyata, N., Design, synthesis, and biological activity of a novel series of human sirtuin-2-selective inhibitors. J Med Chem 2012, 55 (12), 5760-73. [55] Hoffmann, G.; Breitenbucher, F.; Schuler, M.; Ehrenhofer-Murray, A. E., A novel sirtuin 2 (SIRT2) inhibitor with p53-dependent pro-apoptotic activity in non-small cell lung cancer. J Biol Chem 2014, 289 (8), 5208-16.

ACCEPTED MANUSCRIPT

AC C

EP

TE D

M AN U

SC

RI PT

[56] http://www.ebi.ac.uk/. [57] Molecular Operating Environment (MOE); 2012.10, Chemical Computing Group Inc.: Montreal, Canada. 2012. [58] www.rcsb.org. [59] Berman, H. M.; Westbrook, J.; Feng, Z.; Gilliland, G.; Bhat, T. N.; Weissig, H.; Shindyalov, I. N.; Bourne, P. E., The Protein Data Bank. Nucleic Acids Research 2000, 28 (1), 235-242. [60] Wang, J.; Cieplak, P.; Kollman, P. A., How well does a restrained electrostatic potential (RESP) model perform in calculating conformational energies of organic and biological molecules? J Comput Chem 2000, 21 (12), 1049-1074. [61] Avalos, J. L.; Bever, K. M.; Wolberger, C., Mechanism of sirtuin inhibition by nicotinamide: altering the NAD(+) cosubstrate specificity of a Sir2 enzyme. Mol Cell 2005, 17 (6), 855-68. [62] Case, D. A. D., T. A.; Cheatham, III, T. E.; Simmerling, C.; L.; Wang, J. D., R. E.; Luo, R.; Walker, R. C.; Zhang, W.; Merz, K.; M.; et al.; et al. AMBER 12; University of California, S. F., 2012. [63] Glide, Version 5.8; Schrödinger, LLC: New York, 2012. [64] Friesner, R. A.; Banks, J. L.; Murphy, R. B.; Halgren, T. A.; Klicic, J. J.; Mainz, D. T.; Repasky, M. P.; Knoll, E. H.; Shelley, M.; Perry, J. K.; Shaw, D. E.; Francis, P.; Shenkin, P. S., Glide:  A New Approach for Rapid, Accurate Docking and Scoring. 1. Method and Assessment of Docking Accuracy. J Med Chem 2004, 47 (7), 1739-1749. [65] Friesner, R. A.; Murphy, R. B.; Repasky, M. P.; Frye, L. L.; Greenwood, J. R.; Halgren, T. A.; Sanschagrin, P. C.; Mainz, D. T., Extra Precision Glide:  Docking and Scoring Incorporating a Model of Hydrophobic Enclosure for Protein−Ligand Complexes. J Med Chem 2006, 49 (21), 6177-6196. [66] Halgren, T. A.; Murphy, R. B.; Friesner, R. A.; Beard, H. S.; Frye, L. L.; Pollard, W. T.; Banks, J. L., Glide:  A New Approach for Rapid, Accurate Docking and Scoring. 2. Enrichment Factors in Database Screening. J Med Chem 2004, 47 (7), 1750-1759. [67] Duan, Y.; Wu, C.; Chowdhury, S.; Lee, M. C.; Xiong, G.; Zhang, W.; Yang, R.; Cieplak, P.; Luo, R.; Lee, T.; Caldwell, J.; Wang, J.; Kollman, P., A point-charge force field for molecular mechanics simulations of proteins based on condensed-phase quantum mechanical calculations. J Comput Chem 2003, 24 (16), 1999-2012. [68] Lee, M. C.; Duan, Y., Distinguish protein decoys by using a scoring function based on a new AMBER force field, short molecular dynamics simulations, and the generalized born solvent model. Proteins 2004, 55 (3), 620-34. [69] Jakalian, A.; Bush, B. L.; Jack, D. B.; Bayly, C. I., Fast, efficient generation of high-quality atomic charges. AM1-BCC model: I. Method. J Comput Chem 2000, 21 (2), 132-146. [70] Jakalian, A.; Jack, D. B.; Bayly, C. I., Fast, efficient generation of high-quality atomic charges. AM1-BCC model: II. Parameterization and validation. J Comput Chem 2002, 23 (16), 1623-41. [71] Wang, J.; Wolf, R. M.; Caldwell, J. W.; Kollman, P. A.; Case, D. A., Development and testing of a general amber force field. J Comput Chem 2004, 25 (9), 1157-74. [72] http://mayoresearch.mayo.edu/mayo/research/camdl/zinc_protein.cfm. [73] Jorgensen, W. L.; Chandrasekhar, J.; Madura, J. D.; Impey, R. W.; Klein, M. L., Comparison of simple potential functions for simulating liquid water. The Journal of chemical physics 1983, 79 (2), 926-935. [74] Pastor, R. W.; Brooks, B. R.; Szabo, A., An analysis of the accuracy of Langevin and molecular dynamics algorithms. Molecular Physics 1988, 65 (6), 1409-1419. [75] Darden, T.; York, D.; Pedersen, L., Particle mesh Ewald: An N⋅log(N) method for Ewald sums in large systems. The Journal of chemical physics 1993, 98 (12), 10089-10092. [76] Ryckaert, J.-P.; Ciccotti, G.; Berendsen, H. J. C., Numerical integration of the cartesian equations of motion of a system with constraints: molecular dynamics of n-alkanes. Journal of Computational Physics 1977, 23 (3), 327-341. [77] Hawkins, G. D.; Cramer, C. J.; Truhlar, D. G., Pairwise Solute Descreening of Solute

ACCEPTED MANUSCRIPT

AC C

EP

TE D

M AN U

SC

RI PT

Charges from a Dielectric Medium. Chem Phys Lett 1995, 246 (1-2), 122-129. [78] Hawkins, G. D.; Cramer, C. J.; Truhlar, D. G., Parametrized models of aqueous free energies of solvation based on pairwise descreening of solute atomic charges from a dielectric medium. J Phys Chem-Us 1996, 100 (51), 19824-19839. [79] Tsui, V.; Case, D. A., Theory and applications of the generalized Born solvation model in macromolecular Simulations. Biopolymers 2001, 56 (4), 275-291. [80] Onufriev, A.; Bashford, D.; Case, D. A., Exploring protein native states and large-scale conformational changes with a modified generalized born model. Proteins-Structure Function and Bioinformatics 2004, 55 (2), 383-394. [81] Mongan, J.; Simmerling, C.; McCammon, J. A.; Case, D. A.; Onufriev, A., Generalized Born model with a simple, robust molecular volume correction. J Chem Theory Comput 2007, 3 (1), 156169. [82] Wereszczynski, J.; McCammon, J. A., Statistical mechanics and molecular dynamics in evaluating thermodynamic properties of biomolecular recognition. Q Rev Biophys 2012, 45 (1), 125. [83] Hou, T.; Wang, J.; Li, Y.; Wang, W., Assessing the performance of the MM/PBSA and MM/GBSA methods. 1. The accuracy of binding free energy calculations based on molecular dynamics simulations. J Chem Inf Model 2011, 51 (1), 69-82. [84] Wang, R. X.; Lai, L. H.; Wang, S. M., Further development and validation of empirical scoring functions for structure-based binding affinity prediction. J Comput Aid Mol Des 2002, 16 (1), 11-26. [85] Giordanetto, F.; Cotesta, S.; Catana, C.; Trosset, J. Y.; Vulpetti, A.; Stouten, P. F. W.; Kroemer, R. T., Novel scoring functions comprising QXP, SASA, and protein side-chain entropy terms. J Chem Inf Comp Sci 2004, 44 (3), 882-893. [86] Raha, K.; Merz, K. M., A quantum mechanics-based scoring function: Study of zinc ionmediated ligand binding. J Am Chem Soc 2004, 126 (4), 1020-1021. [87] Thompson, D. C.; Humblet, C.; Joseph-McCarthy, D., Investigation of MM-PBSA rescoring of docking poses. J Chem Inf Model 2008, 48 (5), 1081-91. [88] Hayik, S. A.; Dunbrack, R.; Merz, K. M., Mixed Quantum Mechanics/Molecular Mechanics Scoring Function To Predict Protein−Ligand Binding Affinity. J Chem Theory Comput 2010, 6 (10), 3079-3091. [89] Wold, S. E., K.; Geladi, P. , Principal Component Analysis. . Chemom. Intell. Lab. Syst. 1987, 2, 37-52. [90] Charlier, L.; Nespoulous, C.; Fiorucci, S.; Antonczak, S.; Golebiowski, J., Binding free energy prediction in strongly hydrophobic biomolecular systems. Phys Chem Chem Phys 2007, 9 (43), 5761-71. [91] Homeyer, N.; Gohlke, H., FEW: a workflow tool for free energy calculations of ligand binding. J Comput Chem 2013, 34 (11), 965-73. [92] Wang, J.; Morin, P.; Wang, W.; Kollman, P. A., Use of MM-PBSA in reproducing the binding free energies to HIV-1 RT of TIBO derivatives and predicting the binding mode to HIV-1 RT of efavirenz by docking and MM-PBSA. J Am Chem Soc 2001, 123 (22), 5221-30. [93] Hou, T.; Wang, J.; Li, Y.; Wang, W., Assessing the performance of the molecular mechanics/Poisson Boltzmann surface area and molecular mechanics/generalized Born surface area methods. II. The accuracy of ranking poses generated from docking. J Comput Chem 2011, 32 (5), 866-77. [94] Tatiana I. Netzeva, A. P. W., Tom Aldenberg, Romualdo Benigni, Mark T.D.; Cronin, P. G., Joanna S. Jaworska, Scott Kahn, Gilles Klopman, Carol A.; Marchant, G. M., Nina NikolovaJeliazkova, Grace Y. Patlewicz, Roger Perkins,; David W. Roberts, T. W. S., David T. Stanton, Johannes J.M. van de Sandt,; Weida Tong, G. V. a. C. Y., Current Status of Methods for Defining the Applicability Domain of (Quantitative) Structure–Activity Relationships. ATLA 2005, 33, 155173. [95] Stefan Brandmaier, S. N., Iurii Sushko and Igor V. Tetko, From Descriptors to Predicted

ACCEPTED MANUSCRIPT

AC C

EP

TE D

M AN U

SC

RI PT

Properties: Experimental Design by Using Applicability Domain Estimation. ATLA 2013, 41, 33-47. [96] Jaworska, N. N. a. J., Approaches to Measure Chemical Similarity ± a Review. QSAR Comb. Sci. 2003, 22 (9-10), 1006-1026.

ACCEPTED MANUSCRIPT

List of Figure Captions Figure 1. Crystal structure of Sirt3 (white ribbon) complex with 11c (green). The zinc ion is shown as light purple sphere (a). Comparison of the docking solutions of the inhibitors 11c (b), 28 (c) and

RI PT

31 (d) with the location in the corresponding Sirt3 (white ribbon) crystal structures. Only interacting amino acids are shown for clarity. Docking poses are shown in orange whereas the cocrystallized inhibitors are shown in green. Hydrogen bonds are shown as dashed lines in the same

SC

color as the molecule they belong to. Interactions between bridging water molecules and Sirt3 are shown as yellow dashed lines.

M AN U

Figure 2. Superimposition of Sirt3-31 (white sticks) and Sirt1-NAD+-Ex527 (pink sticks). Only amino acids interacting with the inhibitor 31 (shown as green balls and sticks) are depicted. Cofactor NAD+ and the Sirt1 inhibitor Ex527 are not shown for clarity reason. Residues in the catalytic binding site of Sirt3 and Sirt1 structures share a similar conformation except the conserved

TE D

residue Phe157 and Phe273 (cyan sticks) of Sirt3 and Sirt1, respectively. Figure 3. Superimposition of Sirt3-31 (white sticks) and Sirt2-ADPR (pink sticks). Only the amino acids interacting with the inhibitor 31 (shown as green balls and sticks) are depicted. Co-factor analog ADPR is not shown for clarity. Residues in the catalytic binding site of Sirt3 and Sirt2

EP

structures share a similar conformation except the conserved residues Phe157 and Phe96 (cyan

AC C

sticks) of Sirt3 and Sirt2, respectively. Figure 4. Sirt1 (white sticks) in complex with 11c (orange balls and sticks, docking solution). Only the amino acids interacting with the inhibitor 11c are shown for clarity, (a). The molecular surface of the binding pocket is displayed and colored according to the hydrophobicity (green = hydrophobic, magenta = hydrophilic), (b). Figure 5. Sirt2 (white sticks) in complex with 31 (orange balls and sticks, docking solution). Only the amino acids interacting with inhibitor 31 are shown for clarity, (a). The molecular surface of the binding pocket is displayed and colored according to the hydrophobicity (green = hydrophobic,

ACCEPTED MANUSCRIPT magenta = hydrophilic), (b). Figure 6. Docking poses of all inhibitors obtained for Sirt1 (A), Sirt2 (B) and Sirt3 (C). Correlation between GLIDE-XP scores and pIC50 values (using compounds only with exact IC50 data) is shown

RI PT

on the right side for the corresponding Sirtuin. Figure 7. Box plots showing GLIDE-XP score distribution obtained for active and inactive compounds for Sirt1 (A), Sirt2 (B) and Sirt3 (C), respectively.

SC

Figure 8. Correlation between experimental pIC50 values and ∆H values derived from the MMGBSA calculations for Sirt3. ∆H values (given in kcal/mol) were obtained using single

M AN U

conformation after energy minimization in explicit water (A), using 10 conformations obtained from 100 ps MD run (B), single conformation obtained after100 ps equilibration MD run (C). Figure 9. Box plots showing the ∆H score distribution obtained for active and inactive compounds; Sirt1 (A), Sirt2 (B) and Sirt3 (C).

TE D

Figure 10. Superimposition of X-ray structure of Sirt3-11c (Sirt3; white sticks, 11c; green sticks) and Sirt3-31 (pink sticks) in complex with top-ranked docking pose of 11c with GLIDE-XP score (11c; dark pink sticks) and with ∆H score (11c; orange sticks). Only amino acids interacting with

EP

the inhibitor 11c are depicted. Compound 31 is not shown for clarity. Figure 11. Sirt3-31 structure (Sirt3; white sticks, co-crystallized ligand 31; green sticks) in com-

AC C

plex with top-ranked docking pose of 31 with GLIDE-XP score (31; dark pink sticks) and with ∆H score (31; orange sticks). Only amino acids interacting with the inhibitor 31 are depicted.

ACCEPTED MANUSCRIPT Charts Chart 1. Molecular structures of Sirtuin inhibitors used for the current study

RI PT

R2 R1

W R3

( )

nZ

Y

N

W

X

Y

Z

R1

R2

n

11a

S

N

N

CH

H

-(C=O)NH2

2

11b

S

N

N

N

H

-(C=O)NH2

2

H3C

S

N

N

CH

H

-(C=O)NH2

2

11d

S

N

N

N

H

-(C=O)NH2

2

17

S

N

N

S

N

N

CH

N

CH

AC C

19

N

15a

S

N

H

H

EP

S

18

TE D

11c

CH

N

H

-(C=O)NH2

-(C=O)NH2

H3C

NH

NH O

O

NH

NH S O

O

O

NH

2

S O

O

HO

NH

2

S O

O NH

-(C=O)NH2

2

S O

O

CH

H

-(C=O)NH2

2 O

S

N

N

CH

H

-(C=O)NH2

1

34

S

N

N

CH

H

-(C=O)NH2

2

35

S

N

N

N

H

2

20

S

N

N

CH

H

-(C=O)NH2 -(C=O)NH2

23

S

N

N

CH

H

-(C=O)NH2

1

24

S

N

N

N

H

-(C=O)NH2

2

33

R3

M AN U

Cmpd*

SC

X

NH

N

2 -NH(C=O)Me

ACCEPTED MANUSCRIPT 25

S

N

N

CH

H

-(C=O)NH2

2

26

S

N

N

N

H

-(C=O)NH2

27

S

N

N

CH

H

-(C=O)NH2

2 1

28

S

N

N

CH

H

-(C=O)NH2

2

29

S

N

N

CH

H

-(C=O)NH2

3

30

S

N

N

N

H

-(C=O)NH2

2

36

O

N

N

CH

H

-(C=O)NH2

2

37

S

CH

N

CH

H

-(C=O)NH2

2

38

S

N

CH

CH

-(C=O)NH2

2

39

S

N

N

CH

H Me

-(C=O)NH2

40

S

N

N

CH

H

-(C=O)OH

2 2

41

S

N

N

CH

H

-(C=O)NHMe

2

42

S

N

N

CH

H

-H

2

31

S

N

N

CH

H

-(C=O)NH2

2

32

S

N

N

N

H

2

16a

S

N

N

CH

H

-(C=O)NH2 -(C=O)NH2

2

-NH2

21

S

N

N

CH

H

-(C=O)NH2

0

-

RI PT

-NH(C=S)Me

M AN U

SC

-NH(C=O)tBu

-NHSO2Me

TE D

R2

R1

W

R4

Y

X

S

N

22 14

Y

R1

R2

R4

N

H

-(C=O)NH2

EtNH-

H

-(C=O)NH2

Cl-

AC C

Cmpd* W

EP

X

S

N

N

*Cmpd; Compound name

ACCEPTED MANUSCRIPT Tables Table 1. RMSD results (Å) retrieved from cross-docking studies of the inhibitors 11c, 28 and 31 using the Sirt3 crystal structures.

Sirt3-31

(4JSR)

(4JT8)

(4JT9)

11c

1.18

4.35

1.78

28

1.01

0.52

0.72

31

1.03

1.07

1.02

SC

Sirt3-28

AC C

EP

TE D

M AN U

Compound Sirt3-11c

RI PT

Crsytal Structure (PDB ID)

ACCEPTED MANUSCRIPT Table 2. Docking scores, MM-GBSA binding free energies and in vitro data for Sirt1 Activity

M AN U

TE D

EP

AC C

MM-GBSA Score -74.22 -58.88 -73.40 -60.79 -67.55 -64.33 -67.03 -72.77 -65.31 -59.01 -58.21 -62.77 -68.01 -57.85 -58.85 -63.69 -64.13 -57.35 -63.25 -56.34 -60.16 -52.55 -50.15 -55.40 -56.85 -42.89 -54.21 -52.61 -39.83 -57.41 -31.55 -

RI PT

NRot 10 6 10 7 10 8 7 8 10 8 5 6 7 10 5 6 7 7 6 7 6 7 4 4 6 6 2 5 7 3 8 1 6

SC

GLIDE-XP ∆H Compound IC50(µM) pIC50 Score Score 11c 0.0036 8.44 -9.56 -64.22 31 0.0043 8.37 -7.54 -52.88 11a 0.0053 8.28 -8.92 -63.40 19 0.0058 8.24 -7.78 -53.79 17 0.0067 8.17 * 18 0.014 7.85 -9.40 -59.55 28 0.015 7.82 -8.08 -57.33 15a 0.017 7.77 -7.60 -59.03 11d 0.031 7.51 -8.54 -62.77 29 0.042 7.38 -7.85 -57.31 35 0.053 7.28 -8.60 -54.01 25 0.056 7.25 -8.47 -52.21 37 0.059 7.23 -8.27 -55.77 11b 0.067 7.17 -8.29 -58.01 34 0.11 6.96 -8.69 -52.85 20 0.11 6.96 -7.30 -52.85 30 0.12 6.92 -8.76 -56.69 39 0.18 6.74 -7.70 -57.13 32 0.38 6.42 -7.25 -51.35 36 0.47 6.33 -7.93 -56.25 26 0.49 6.31 -7.27 -50.34 38 0.49 6.31 -6.07 -53.16 33 0.73 6.14 -7.38 -48.55 16a 1.6 5.80 -8.70 -46.15 24 1.6 5.80 -7.14 -49.40 27 4.3 5.37 -7.61 -50.85 21 15 4.82 -7.55 -40.89 23 16 4.80 -7.16 -49.21 40 >50 -7.98 -45.61 22 >50 -6.84 -36.83 41 >50 -6.29 -49.41 14 >50 -5.79 -30.55 42 >50 * *GLIDE docking did not return back any pose for the compound.

ACCEPTED MANUSCRIPT Table 3. Docking scores, MM-GBSA binding free energies and in vitro data for Sirt2 Activity MM-GBSA Score -68.88 -79.07 -79.47 -73.19 -76.34 -79.43 -72.71 -69.67 -69.97 -78.40 -74.68 -60.89 -66.38 -69.17 -62.30 -65.89 -66.11 -63.81 -57.01 -71.15 -65.73 -63.25 -65.13 -63.02 -66.79 -61.57 -49.30 -61.20 -67.57 -46.59 -68.43 -49.60 -38.96

RI PT

NRot 6 10 10 10 8 10 8 7 7 10 8 7 6 7 5 6 6 5 4 7 6 7 6 4 7 6 2 5 8 3 7 6 1

SC

∆H Score -62.88 -69.07 -69.47 -63.19 -68.34 -69.43 -64.71 -62.67 -62.97 -68.40 -66.68 -53.89 -60.38 -62.17 -57.30 -59.89 -60.11 -58.81 -53.01 -64.15 -59.73 -56.25 -59.13 -59.02 -59.79 -55.57 -47.30 -56.20 -59.57 -43.59 -61.43 -43.60 -37.96

M AN U

GLIDE-XP Score -8.64 -9.14 -9.16 -9.82 -10.59 -10.48 -10.29 -10.76 -10.51 -10.07 -10.24 -8.52 -8.38 -9.03 -9.23 -10.17 -9.83 -9.53 -9.59 -9.92 -9.55 -9.58 -10.10 -9.57 -9.22 -9.55 -8.95 -9.63 -8.29 -7.94 -10.69 -7.54 -7.54

TE D

pIC50 8.96 8.89 8.74 8.57 8.36 8.32 8.27 8.19 8.00 8.00 7.80 7.77 7.64 7.55 7.54 7.52 7.28 7.17 6.96 6.89 6.89 6.82 6.77 6.66 6.52 6.13 5.85 5.55 5.00 4.31 -

EP

IC50(µM) 0.0011 0.0013 0.0018 0.0027 0.0044 0.0048 0.0054 0.0065 0.010 0.010 0.016 0.017 0.023 0.028 0.029 0.030 0.053 0.067 0.11 0.13 0.13 0.15 0.17 0.22 0.30 0.74 1.4 2.8 10 49 >50 >50 >50

AC C

Compound 31 11a 17 11c 18 11d 15a 19 28 11b 29 30 20 37 35 25 32 34 16a 39 26 36 24 33 38 27 21 23 41 22 40 42 14

ACCEPTED MANUSCRIPT Table 4. Docking scores, MM-GBSA binding free energies and in vitro data for Sirt3 Activity MM-GBSA Score -74.42 -65.82 -75.27 -60.19 -65.41 -72.30 -75.07 -73.14 -61.36 -67.47 -60.10 -66.17 -60.64 -59.35 -57.07 -68.09 -65.53 -65.82 -66.07 -58.98 -50.90 -60.14 -59.21 -59.20 -63.85 -63.36 -58.04 -43.57 -62.80 -58.00 -49.98 -39.01 -34.25

RI PT

NRot 10 10 10 6 7 8 10 10 5 8 5 7 6 6 4 8 7 7 7 6 4 6 6 6 7 7 5 2 8 7 6 3 1

SC

∆H Score -64.42 -55.82 -65.27 -54.19 -58.41 -64.30 -65.07 -63.14 -56.36 -59.47 -55.10 -59.17 -54.64 -53.35 -53.07 -60.09 -58.53 -58.82 -59.07 -52.98 -46.90 -54.14 -53.21 -53.20 -56.85 -56.36 -53.04 -41.57 -54.80 -51.00 -43.98 -36.01 -33.25

M AN U

GLIDE-XP Score -12.21 -10.24 -11.93 -11.63 -11.36 -12.35 -10.49 -11.45 -10.40 -10.04 -10.38 -12.07 -11.87 -11.44 -10.22 -11.86 -11.53 -11.92 -9.04 -11.33 -11.27 -11.10 -11.23 -10.70 -11.42 -9.81 -10.85 -9.23 -8.72 -10.18 -8.62 -8.31 -8.14

TE D

pIC50 8.49 8.40 8.31 8.14 8.12 7.89 7.64 7.62 7.55 7.54 7.52 7.48 7.41 7.25 7.21 7.17 7.04 6.96 6.82 6.55 6.52 6.43 6.19 6.11 5.89 5.68 5.36 5.07 -

EP

IC50(µM) 0.0032 0.0040 0.0049 0.0072 0.0076 0.013 0.023 0.024 0.028 0.029 0.030 0.033 0.039 0.056 0.061 0.067 0.092 0.11 0.15 0.28 0.30 0.37 0.65 0.77 1.3 2.1 4.4 8.6 >50 >50 >50 >50 >50

AC C

Compound 17 11c 11a 31 19 18 11d 11b 34 15a 35 28 25 20 33 29 30 37 39 26 16a 32 27 24 36 38 23 21 41 40 42 22 14

ACCEPTED MANUSCRIPT Table 5. Summary of the statistics obtained for MM-PBSA and MM-GBSA models for Sirt3 md 0.21 (n=28) 0.80 31 = 2.02 36 = 2.02 38 = 2.17 11c = 2.22 0.50 (n=24) 0.59 0.43

Method

Statistics R2 RMSE Outliers (XZ-score)

min 0.25 (n=28) 0.79 11c = 2.15 23 = 2.59

R2 RMSE q2LOO

0.36 (n=26) 0.67 0.24

md 0.42 (n=28) 0.69 31 = 2.13 23 = 2.14 38 = 2.19 11c = 2.26 0.60 (n=24) 0.50 0.54

GB-1

GB-1*

Method

Statistics R2 RMSE Outliers (XZ-score)

min 0.22 (n=28) 0.80 23 = 2.17 31 = 2.26 11c = 2.31 0.42 (n=25) 0.63 0.31

md 0.30 (n=28) 0.76 23 = 2.33 38 = 2.48

GB-2

min 0.25 (n=28) 0.79 38 = 2.05 23 = 2.06 11c = 2.10 0.34 (n=25) 0.65 0.21

md 0.30 (n=28) 0.75 11c = 2.06 23 = 2.28 38 = 2.36 0.45 (n=25) 0.60 0.35

AC C

GB-2*

R2 RMSE q2LOO

TE D

PB-1*

EP

PB-1

R RMSE Outliers (XZ-score)

Method

GB-5

GB-5*

Statistics R2 RMSE Outliers (XZ-score) R2 RMSE q2LOO

md-final 0.25 (n=28) 0.79 31 = 2.18 38 = 2.30 17 = 2.40

RI PT

R2 RMSE q2LOO

min 0.20 (n=28) 0.81 31 = 2.04 36 = 2.10 38 = 2.17 11c = 2.29 0.53 (n=24) 0.57 0.45

SC

Statistics 2

M AN U

Method

0.40 (n=26) 0.63 0.30

0.49 (n=25) 0.60 0.41

md-final 0.44 (n=28) 0.68 31 = 2.18 38 = 2.31

0.54 (n=26) 0.59 0.48

md-final 0.39 (n=28) 0.71 31 = 2.09 38 = 2.20 11c = 2.21 0.55 (n=25) 0.57 0.48

md-final 0.37 (n=28) 0.72 23 = 2.22 38 = 2.64 0.47 (n=26) 0.60 0.38

ACCEPTED MANUSCRIPT Table 5. Summary of the statistics obtained for MM-PBSA and MM-GBSA models for Sirt3 (continued) min md md-final 0.34 (n=28) 0.34 (n=28) 0.34 (n=28) 0.74 0.74 0.73 36 = 2.04 36 = 2.04 38 = 2.68 GB-8 23 = 2.07 23 = 2.11 11c = 2.23 11c = 2.19 38 = 2.46 38 = 2.53 R2 0.57 (n=24) 0.57 (n=24) 0.42 (n=27) RMSE 0.51 0.51 0.67 GB-8* q2LOO 0.47 0.48 0.34 * Statistical values calculated for the MM-GBSA method after the omission of compounds

RI PT

Statistics R2 RMSE Outliers (XZ-score)

SC

Method

identified as outliers based on their XZ-score. n; number of compounds used to analyse the

AC C

EP

TE D

M AN U

statistics.

ACCEPTED MANUSCRIPT Table 6. Summary of the statistics obtained for MM-GBSA models for Sirt1 Statistics

R2 RMSE q2LOO

min 0.47 (n=27) 0.73 19 = 2.00 23 = 2.17 31 = 3.16 0.63 (n=24) 0.54 0.56

md 0.60 (n=27) 0.63 23 = 2.48 19 = 2.54 31 = 3.19 0.79 (n=24) 0.41 0.76

Statistics R2 RMSE Outliers (XZ-score)

min 0.44 (n=27) 0.75 23 = 2.71 31 = 4.13

md 0.44 (n=27) 0.75 23 = 2.72 31 = 4.44

2

GB-1

GB-1*

SC

md-final 0.39 (n=27) 0.78 GB-8 19 = 2.01 23 = 2.73 31 = 4.27 R2 0.67 (n=25) 0.69 (n=25) 0.74 (n=24) RMSE 0.52 0.51 0.45 GB-8* q2LOO 0.61 0.64 0.69 * Statistical values calculated for the MM-GBSA method after the omission of compounds

M AN U

Method

R RMSE Outliers (XZ-score)

md-final 0.58 (n=27) 0.65 23 = 2.32 19 = 2.49 31 = 3.12 0.75 (n=24) 0.44 0.71

RI PT

Method

AC C

EP

statistics.

TE D

identified as outliers based on their XZ-score. n; number of compounds used to analyse the

ACCEPTED MANUSCRIPT Table 7. Summary of the statistics obtained for MM-GBSA models for Sirt2 Statistics 2

GB-1

GB-1*

R2 RMSE q2LOO

md-final 0.61 (n=30) 0.69 28 = 2.06 30 = 3.32 0.71 (n=28) 0.57 0.67

min md md-final 0.46 (n=30) 0.54 (n=30) 0.55 (n=30) 0.82 0.76 0.75 GB-8 30 = 2.21 31 = 2.08 30 = 2.28 31 = 2.73 30 = 2.36 41 = 3.40 41 = 3.09 41 = 3.36 R2 0.69 (n=27) 0.72 (n=27) 0.68 (n=28) RMSE 0.58 0.55 0.60 GB-8* q2LOO 0.62 0.66 0.63 * Statistical values calculated for the MM-GBSA method after the omission of compounds

SC

Statistics R2 RMSE Outliers (XZ-score)

md 0.59 (n=30) 0.72 31 = 2.00 30 = 2.24 41 = 3.74 0.76 (n=27) 0.51 0.72

M AN U

Method

R RMSE Outliers (XZ-score)

min 0.54 (n=30) 0.76 23 = 2.26 35 = 2.39 22 = 2.98 0.63 (n=27) 0.58 0.58

RI PT

Method

AC C

EP

statistics.

TE D

identified as outliers based on their XZ-score. n; number of compounds used to analyse the

ACCEPTED MANUSCRIPT Table 8. ∆H Scores obtained for the outliers using the top-ranked poses with GLIDE-XP scoring and MM-GBSA rescoring. ∆H Score (kcal/mol) A B -53.79 -58.84 19 Sirt1 -52.88 -55.73 31 -53.89 -57.84 30 Sirt2 -62.88 -62.88 31* -55.82 -60.23 11c Sirt3 -54.11 -56.68 31 *GLIDE docking returned only one pose for compound 31. A; ∆H score obtained for the top-ranked

RI PT

Outlier

AC C

EP

TE D

M AN U

SC

docking pose, B; ∆H score obtained for the lowest binding energy pose.

ACCEPTED MANUSCRIPT Table 9. Summary of the statistics obtained for the MM-GBSA scores versus pIC50 correlation plots for the refined MMGBSA Method Sirt1 Sirt2 Sirt3 0.58 (n=27) 0.57 (n=30) 0.43 (n=28) MM-GBSA 0.65 0.73 0.68 (entropy 19 = 2.32 31 = 2.23 23 = 2.08 included) 23 = 2.44 41 = 4.07 31 = 2.16 31 = 3.09 38 = 2.27 2 R 0.75 (n=24) 0.72 (n=28) 0.54 (n=25) MM-GBSA* RMSE 0.44 0.54 0.55 (entropy included) q2LOO 0.71 0.66 0.47 * Statistical values calculated for the MMGBSA method after omitting compounds stated as outliers

RI PT

Statistics R2 RMSE Outliers (XZ-score)

SC

Method

AC C

EP

TE D

M AN U

based on their XZ-score. n; number of compounds used to analyse the statistics.

ACCEPTED MANUSCRIPT

EP

TE D

M AN U

SC

RI PT

Figures

Figure 1. Crystal structure of Sirt3 (white ribbon) complex with 11c (green). The zinc ion is shown

AC C

as light purple sphere (a). Comparison of the docking solutions of the inhibitors 11c (b), 28 (c) and 31 (d) with the location in the corresponding Sirt3 (white ribbon) crystal structures. Only interacting amino acids are shown for clarity. Docking poses are shown in orange whereas the cocrystallized inhibitors are shown in green. Hydrogen bonds are shown as dashed lines in the same color as the molecule they belong to. Interactions between bridging water molecules and Sirt3 are shown as yellow dashed lines.

SC

RI PT

ACCEPTED MANUSCRIPT

M AN U

Figure 2. Superimposition of Sirt3-31 (white sticks) and Sirt1-NAD+-Ex527 (pink sticks). Only amino acids interacting with the inhibitor 31 (shown as green balls and sticks) are depicted. Cofactor NAD+ and the Sirt1 inhibitor Ex527 are not shown for clarity reason. Residues in the catalytic binding site of Sirt3 and Sirt1 structures share a similar conformation except the conserved

AC C

EP

TE D

residue Phe157 and Phe273 (cyan sticks) of Sirt3 and Sirt1, respectively.

SC

RI PT

ACCEPTED MANUSCRIPT

M AN U

Figure 3. Superimposition of Sirt3-31 (white sticks) and Sirt2-ADPR (pink sticks). Only the amino acids interacting with the inhibitor 31 (shown as green balls and sticks) are depicted. Co-factor analog ADPR is not shown for clarity. Residues in the catalytic binding site of Sirt3 and Sirt2 structures share a similar conformation except the conserved residues Phe157 and Phe96 (cyan

AC C

EP

TE D

sticks) of Sirt3 and Sirt2, respectively.

RI PT

ACCEPTED MANUSCRIPT

Figure 4. Sirt1 (white sticks) in complex with 11c (orange balls and sticks, docking solution). Only

SC

the amino acids interacting with the inhibitor 11c are shown for clarity, (a). The molecular surface

AC C

EP

TE D

hydrophobic, magenta = hydrophilic), (b).

M AN U

of the binding pocket is displayed and colored according to the hydrophobicity (green =

RI PT

ACCEPTED MANUSCRIPT

Figure 5. Sirt2 (white sticks) in complex with 31 (orange balls and sticks, docking solution). Only

SC

the amino acids interacting with inhibitor 31 are shown for clarity, (a). The molecular surface of the

AC C

EP

TE D

magenta = hydrophilic), (b).

M AN U

binding pocket is displayed and colored according to the hydrophobicity (green = hydrophobic,

ACCEPTED MANUSCRIPT -10

Phe297

R² = 0.26

Ile270

W522 Ile347

-9 -8 -7 -6

Asp348 Ala262

His363

(A)

GLIDE-XP Score

Phe273

-5 3

RI PT

Phe414

4

5

6

7

8

9 10

-11

Phe119

R² = 0.18

Ile93

His187

W2024

Ile169 Asp170

(B)

Ala85

-9

-8

-7

3

7

8

9 10

Ile154

W522 Ile230 Ala146

8

9 10

-12

GLIDE-XP Score

EP AC C

6

R² = 0.16

Phe157

Asp231

5

-13

Phe180

His248

4

pIC50

TE D

Phe294

(C)

-10

M AN U

Phe96

GLIDE-XP Score

Phe235

SC

pIC50

-11 -10 -9 -8 3

4

5

6

7

pIC50

Figure 6. Docking poses of all inhibitors obtained for Sirt1 (A), Sirt2 (B) and Sirt3 (C). Correlation

between GLIDE-XP scores and pIC50 values (using compounds only with exact IC50 data) is shown on the right side for the corresponding Sirtuin.

(A)

(B)

RI PT

ACCEPTED MANUSCRIPT

(C)

Figure 7. Box plots showing GLIDE-XP score distribution obtained for active and inactive

AC C

EP

TE D

M AN U

SC

compounds for Sirt1 (A), Sirt2 (B) and Sirt3 (C), respectively.

ACCEPTED MANUSCRIPT

(B)

(C)

EP

(A)

TE D

M AN U

GB-8 METHOD

SC

RI PT

GB-1 METHOD

Figure 8. Correlation between experimental pIC50 values and ∆H values derived from the MM-

AC C

GBSA calculations for Sirt3. ∆H values (given in kcal/mol) were obtained using single conformation after energy minimization in explicit water (A), using 10 conformations obtained from 100 ps MD run (B), single conformation obtained after100 ps equilibration MD run (C).

(A)

(B)

RI PT

ACCEPTED MANUSCRIPT

(C)

Figure 9. Box plots showing the ∆H score distribution obtained for active and inactive compounds;

AC C

EP

TE D

M AN U

SC

Sirt1 (A), Sirt2 (B) and Sirt3 (C).

SC

RI PT

ACCEPTED MANUSCRIPT

M AN U

Figure 10. Superimposition of X-ray structure of Sirt3-11c (Sirt3; white sticks, 11c; green sticks)

and Sirt3-31 (pink sticks) in complex with top-ranked docking pose of 11c with GLIDE-XP score (11c; dark pink sticks) and with ∆H score (11c; orange sticks). Only amino acids interacting with

AC C

EP

TE D

the inhibitor 11c are depicted. Compound 31 is not shown for clarity.

SC

RI PT

ACCEPTED MANUSCRIPT

Figure 11. Sirt3-31 structure (Sirt3; white sticks, co-crystallized ligand 31; green sticks) in com-

M AN U

plex with top-ranked docking pose of 31 with GLIDE-XP score (31; dark pink sticks) and with ∆H

AC C

EP

TE D

score (31; orange sticks). Only amino acids interacting with the inhibitor 31 are depicted.

ACCEPTED MANUSCRIPT

HIGHLIGHTS

EP

TE D

M AN U

SC

RI PT

Robust QSAR models were derived for sirtuin inhibitors acting on Sirt1-3. MM-GBSA method showed better agreement with the activity data than docking scores. The developed MM-GBSA method represents a time-efficient post-docking filter. The best overall performance was observed using the igb = 1 solvation model. Promising results were obtained using the fairly new igb = 8 solvation model.

AC C

• • • • •

ACCEPTED MANUSCRIPT Supporting Information for the Manuscript

Institute of Pharmacy, MLU Halle-Wittenberg, 06120 Halle (Saale), Germany

M AN U



SC

Berin Karaman† and Wolfgang Sippl†

RI PT

Docking and binding free energy calculations of sirtuin inhibitors

Table of Content

Page S3. Figure S1. RMSD plot for MD simulation of Sirt3-17.

TE D

Page S4. Figure S2. Superimposition of Sirt3-31 structure with Sirt1 structures. Page S5. Figure S3. RMSD plot for MD simulation of Sirt1-11c.

EP

Page S6. Figure S4. Superimposition of Sirt3-31 structure with Sirt2 structures. Page S7. Figure S5. RMSD plot for MD simulation of Sirt2-31.

AC C

Page S8. Table S1. MM-PBSA and MM-GBSA binding free energies (kcal mol-1) calculated with the final conformation after minimization in explicit water (min) and in vitro data for Sirt3. Page S9. Table S2. MM-PBSA and MM-GBSA binding free energies (kcal mol-1) calculated with 10 conformations obtained in 100 ps equilibration MD (md) and in vitro data for Sirt3. Page S10. Table S3. MM-PBSA and MM-GBSA binding free energies (kcal mol-1) calculated with the final conformation after 100 ps equilibration MD (md-final) and in vitro data for Sirt3.

ACCEPTED MANUSCRIPT Page S11. Table S4. MM-GBSA binding free energies (kcal mol-1) and in vitro data for Sirt1. Page S12. Table S5. MM-GBSA binding free energies (kcal mol-1) and in vitro data for Sirt2. Page S13. Table S6. Summary of the statistics obtained for the MM-GBSA model based on the

RI PT

lowest binding free energy scores.

Page S14. Figure S6. Principal component analysis (PCA) plot obtained for the studied sirtuin

AC C

EP

TE D

M AN U

SC

inhibitors by using simple 1D descriptors.

M AN U

SC

RI PT

ACCEPTED MANUSCRIPT

Figure S1. RMSD plot for MD simulation of Sirt3-17. The black line represents RMSD value

TE D

calculated for backbone heavy atoms of the whole protein. The red line represents the RMSD value calculated for backbone heavy atoms of pocket residues only. The green line represents the RMSD

AC C

EP

value calculated for heavy atoms of the inhibitor, 17.

SC

RI PT

ACCEPTED MANUSCRIPT

M AN U

Figure S2. Superimposition of Sirt3-31 structure with Sirt1 structures. Alpha carbon atoms are considered for RMSD calculations. Bound conformations of Sirt1, Sirt1-ADPR and Sirt1-NAD+Ex527 structures, showed more resemblance to Sirt3-31 structure (RMSD of all alpha carbon

AC C

EP

TE D

atoms; 0.94 Å and 1.03 Å, respectively) than Sirt1-Apo (RMSD of all alpha carbon atoms; 1.54 Å).

M AN U

SC

RI PT

ACCEPTED MANUSCRIPT

TE D

Figure S3. RMSD plot for MD simulation of Sirt1-11c. The black line represents RMSD value calculated for backbone heavy atoms of the whole protein. The red line represents the RMSD value calculated for backbone heavy atoms of pocket residues only. The green line represents the RMSD

AC C

EP

value calculated for heavy atoms of the inhibitor, 11c.

M AN U

SC

RI PT

ACCEPTED MANUSCRIPT

Figure S4. Superimposition of Sirt3-31 structure with Sirt2 structures. Alpha carbon atoms are considered for RMSD calculations. Sirt2-ADPR structure showed more resemblance to Sirt3-31 structure (RMSD of all alpha carbon atoms; 1.06 Å) than Sirt1-Apo (RMSD of all alpha carbon

AC C

EP

TE D

atoms; 1.48 Å).

M AN U

SC

RI PT

ACCEPTED MANUSCRIPT

TE D

Figure S5. RMSD plot for MD simulation of Sirt2-31. The black line represents RMSD value calculated for backbone heavy atoms of the whole protein. The red line represents the RMSD value calculated for backbone heavy atoms of pocket residues only. The green line represents the RMSD

AC C

EP

value calculated for heavy atoms of the inhibitor, 31.

ACCEPTED MANUSCRIPT Table S1. MM-PBSA and MM-GBSA binding free energies (kcal mol-1) calculated with the final conformation after minimization in explicit water (min) and in vitro data for Sirt3

EP

AC C

GB-2 -58.40 -49.22 -59.08 -46.91 -56.85 -61.08 -61.60 -61.22 -55.62 -54.15 -48.13 -57.74 -53.36 -51.26 -47.78 -54.99 -58.45 -58.52 -54.20 -48.52 -41.87 -54.16 -51.80 -47.37 -55.26 -51.69 -51.14 -40.29 -31.76 -34.32 -50.70 -50.55 -39.34

GB-5 -58.40 -52.09 -59.08 -50.62 -56.85 -61.08 -61.60 -61.22 -55.62 -57.69 -48.13 -57.74 -53.36 -51.26 -51.78 -58.41 -58.45 -58.52 -57.65 -48.52 -41.87 -54.16 -51.80 -51.81 -55.26 -55.02 -51.14 -40.29 -31.76 -34.32 -50.70 -50.55 -41.17

GB-8 -53.92 -46.55 -56.88 -45.97 -51.90 -58.96 -53.68 -49.75 -47.89 -48.62 -44.23 -51.74 -48.46 -46.52 -45.90 -53.06 -49.40 -53.28 -49.15 -45.59 -39.32 -39.16 -44.72 -44.47 -49.85 -50.10 -44.59 -34.72 -28.55 -30.09 -44.98 -48.09 -37.30

RI PT

GB-1 -59.79 -49.19 -54.62 -50.88 -52.50 -61.49 -57.79 -63.45 -53.61 -56.19 -50.21 -43.93 -48.51 -48.95 -51.88 -46.81 -48.76 -50.65 -51.80 -52.10 -47.92 -55.94 -45.60 -49.94 -48.35 -45.97 -52.22 -42.79 -34.63 -35.42 -40.22 -37.99 -35.84

SC

PB-1 -43.93 -40.87 -52.61 -39.97 -48.20 -50.92 -52.14 -42.90 -43.11 -43.28 -45.32 -49.34 -45.00 -42.98 -43.19 -49.48 -48.77 -48.61 -44.10 -42.77 -37.63 -42.05 -40.99 -40.88 -48.01 -46.44 -40.02 -31.42 -35.47 -43.91 -33.55 -26.49 -21.71

M AN U

17 11c 11a 31 19 18 11d 11b 34 15a 35 28 25 20 33 29 30 37 39 26 16a 32 27 24 36 38 23 21 14 22 40 41 42

Activity IC50 (µM) pIC50 0.0032 8.49 0.0040 8.40 0.0049 8.31 0.0072 8.14 0.0076 8.12 0.0130 7.89 0.0230 7.64 0.0240 7.62 0.0280 7.55 0.0290 7.54 0.0300 7.52 0.0330 7.48 0.0390 7.41 0.0560 7.25 0.0610 7.21 0.0670 7.17 0.0920 7.04 0.1100 6.96 0.1500 6.82 0.2800 6.55 0.3000 6.52 0.3700 6.43 0.6500 6.19 0.7700 6.11 1.3000 5.89 2.1000 5.68 4.4000 5.36 8.6000 5.07 >50 >50 >50 >50 >50 -

TE D

Compound

ACCEPTED MANUSCRIPT Table S2. MM-PBSA and MM-GBSA binding free energies (kcal mol-1) calculated with 10 conformations obtained in 100 ps equilibration MD (md) and in vitro data for Sirt3

EP

AC C

GB-2 -56.27 -50.45 -56.85 -49.22 -54.98 -57.77 -55.07 -59.44 -49.07 -55.89 -46.80 -52.36 -48.06 -46.28 -49.84 -52.87 -56.67 -52.13 -51.76 -45.79 -43.45 -51.01 -49.21 -45.30 -50.53 -53.44 -49.66 -38.83 -28.52 -31.07 -48.63 -48.19 -40.25

GB-5 -59.96 -50.45 -60.29 -49.22 -54.98 -57.77 -55.07 -59.44 -49.07 -55.89 -46.80 -52.36 -48.06 -46.28 -49.84 -52.87 -56.67 -55.94 -55.34 -45.79 -43.45 -51.01 -49.21 -45.30 -50.53 -53.44 -49.66 -38.83 -28.52 -33.03 -48.63 -51.15 -40.25

GB-8 -51.75 -44.44 -54.66 -43.25 -50.71 -55.53 -51.94 -47.77 -45.49 -45.94 -41.84 -49.51 -46.29 -44.27 -44.03 -50.66 -46.87 -50.56 -46.56 -42.83 -37.49 -35.96 -42.23 -42.49 -47.58 -48.14 -42.60 -32.90 -45.85 -42.75 -36.28 -28.33 -27.71

RI PT

GB-1 -64.42 -55.82 -65.27 -54.19 -58.41 -64.30 -65.07 -63.14 -56.36 -59.47 -55.10 -59.17 -54.64 -53.35 -53.07 -60.09 -58.53 -58.82 -59.07 -52.98 -46.90 -54.14 -53.21 -53.20 -56.85 -56.36 -53.04 -41.57 -33.25 -36.01 -51.00 -54.80 -43.98

SC

PB-1 -41.35 -39.05 -50.09 -37.90 -44.86 -47.23 -50.04 -40.37 -40.58 -40.74 -42.98 -47.05 -42.75 -41.13 -40.96 -46.63 -46.33 -45.78 -42.07 -40.37 -35.32 -39.80 -38.03 -39.42 -44.90 -43.83 -38.08 -29.10 -22.09 -24.78 -33.01 -41.23 -32.66

M AN U

17 11c 11a 31 19 18 11d 11b 34 15a 35 28 25 20 33 29 30 37 39 26 16a 32 27 24 36 38 23 21 14 22 40 41 42

Activity IC50 (µM) pIC50 0.0032 8.49 0.0040 8.40 0.0049 8.31 0.0072 8.14 0.0076 8.12 0.0130 7.89 0.0230 7.64 0.0240 7.62 0.0280 7.55 0.0290 7.54 0.0300 7.52 0.0330 7.48 0.0390 7.41 0.0560 7.25 0.0610 7.21 0.0670 7.17 0.0920 7.04 0.1100 6.96 0.1500 6.82 0.2800 6.55 0.3000 6.52 0.3700 6.43 0.6500 6.19 0.7700 6.11 1.3000 5.89 2.1000 5.68 4.4000 5.36 8.6000 5.07 >50 >50 >50 >50 >50 -

TE D

Compound

ACCEPTED MANUSCRIPT Table S3. MM-PBSA and MM-GBSA binding free energies (kcal mol-1) calculated with the final conformation after 100 ps equilibration MD (md-final) and in vitro data for Sirt3

EP

AC C

GB-2 -58.86 -49.47 -58.01 -47.66 -55.75 -56.89 -53.29 -60.83 -52.84 -52.35 -46.88 -53.24 -46.68 -46.62 -45.48 -56.64 -52.17 -55.70 -53.36 -46.27 -41.25 -51.05 -44.12 -46.58 -49.66 -50.51 -45.37 -35.53 -27.86 -30.11 -47.98 -50.76 -37.77

GB-5 -58.86 -52.43 -58.01 -51.46 -55.75 -61.17 -53.29 -60.83 -52.84 -55.95 -46.88 -53.24 -46.68 -46.62 -49.44 -56.64 -52.17 -55.70 -57.03 -46.27 -44.66 -51.05 -44.12 -46.58 -49.66 -54.19 -48.81 -38.04 -29.69 -32.20 -47.98 -50.76 -39.46

GB-8 -50.53 -46.26 -55.21 -44.34 -51.07 -57.71 -51.06 -48.53 -44.71 -46.75 -40.93 -50.74 -44.37 -42.85 -43.57 -50.69 -49.47 -50.48 -47.94 -43.44 -38.41 -35.78 -43.60 -43.72 -45.85 -49.13 -40.82 -33.02 -27.21 -28.29 -41.90 -45.21 -35.42

RI PT

GB-1 -63.75 -57.67 -65.25 -54.73 -59.07 -63.22 -64.75 -63.97 -58.16 -59.86 -55.00 -60.14 -53.08 -53.52 -54.85 -60.51 -58.32 -60.59 -58.76 -53.71 -47.46 -55.74 -53.58 -54.30 -56.48 -57.02 -52.08 -41.67 -33.02 -35.17 -50.78 -54.06 -43.25

SC

PB-1 -39.95 -42.90 -48.53 -38.08 -45.85 -46.93 -50.22 -41.55 -43.76 -41.21 -43.13 -49.07 -39.50 -42.61 -43.94 -47.05 -45.25 -46.82 -43.36 -40.61 -36.26 -42.38 -37.62 -41.21 -44.62 -44.63 -34.92 -29.25 -22.82 -26.53 -33.09 -39.17 -30.51

M AN U

17 11c 11a 31 19 18 11d 11b 34 15a 35 28 25 20 33 29 30 37 39 26 16a 32 27 24 36 38 23 21 14 22 40 41 42

Activity IC50 (µM) pIC50 0.0032 8.49 0.0040 8.40 0.0049 8.31 0.0072 8.14 0.0076 8.12 0.0130 7.89 0.0230 7.64 0.0240 7.62 0.0280 7.55 0.0290 7.54 0.0300 7.52 0.0330 7.48 0.0390 7.41 0.0560 7.25 0.0610 7.21 0.0670 7.17 0.0920 7.04 0.1100 6.96 0.1500 6.82 0.2800 6.55 0.3000 6.52 0.3700 6.43 0.6500 6.19 0.7700 6.11 1.3000 5.89 2.1000 5.68 4.4000 5.36 8.6000 5.07 >50 >50 >50 >50 >50 -

TE D

Compound

ACCEPTED MANUSCRIPT Table S4. MM-GBSA binding free energies (kcal mol-1) and in vitro data for Sirt1

min -52.78 -37.52 -51.17 -44.57 -52.57 -46.06 -48.41 -51.07 -45.97 -41.08 -42.56 -46.95 -45.17 -42.38 -44.46 -44.24 -48.61 -35.87 -46.53 -41.81 -43.02 -39.78 -37.29 -37.88 -39.20 -32.57 -41.04 -37.66 -29.67 -40.07 -25.93 -

GB-8 md md-final -50.25 -49.27 -35.35 -34.73 -48.56 -47.89 -42.14 -41.04 -49.74 -48.70 -44.35 -46.02 -46.33 -45.14 -48.20 -49.20 -43.81 -44.57 -39.76 -40.42 -41.47 -42.57 -44.71 -46.45 -43.30 -43.51 -40.62 -40.36 -42.68 -43.94 -42.40 -44.06 -45.90 -45.89 -34.08 -34.92 -44.00 -44.67 -39.99 -37.56 -41.65 -42.43 -38.28 -38.57 -35.60 -33.95 -36.99 -36.11 -36.90 -37.48 -30.90 -30.86 -39.30 -39.87 -36.08 -37.49 -27.89 -28.95 -38.06 -38.04 -24.43 -23.60 -

AC C

EP

TE D

M AN U

SC

RI PT

Activity GB-1 IC50 (µM) pIC50 min md md-final 8.44 -50.88 -64.22 -64.64 11c 0.0036 8.37 -43.54 -52.88 -52.45 31 0.0043 8.28 -53.29 -63.40 -62.74 11a 0.0053 8.24 -46.28 -53.79 -53.47 19 0.0058 17* 0.0067 8.17 7.85 -56.29 -59.55 -58.13 18 0.0140 7.82 -47.00 -57.33 -58.36 28 0.0150 7.77 -46.38 -59.03 -57.93 15a 0.0170 7.51 -54.83 -62.77 -63.85 11d 0.0310 7.38 -46.62 -57.31 -57.96 29 0.0420 7.28 -47.76 -54.01 -53.31 35 0.0530 7.25 -45.58 -52.21 -52.81 25 0.0560 7.23 -47.39 -55.77 -54.33 37 0.0590 7.17 -50.87 -58.01 -57.96 11b 0.0670 6.96 -47.25 -52.85 -53.32 34 0.1100 6.96 -44.28 -52.85 -51.86 20 0.1100 6.92 -48.60 -56.69 -58.46 30 0.1200 6.74 -48.16 -57.13 -56.80 39 0.1800 6.42 -44.70 -51.35 -51.74 32 0.3800 6.33 -44.38 -56.25 -56.66 36 0.4700 6.31 -42.59 -50.34 -51.33 26 0.4900 6.31 -41.05 -53.16 -52.88 38 0.4900 6.14 -44.76 -48.55 -47.88 33 0.7300 5.80 -41.00 -46.15 -45.95 16a 1.6000 5.80 -43.76 -49.40 -49.09 24 1.6000 5.37 -44.58 -50.85 -50.80 27 4.3000 4.82 -36.72 -40.89 -39.58 21 15.0000 4.80 -41.86 -49.21 -48.46 23 16.0000 >50 -36.12 -45.61 -45.76 40 >50 -31.35 -36.83 -36.04 22 >50 -33.48 -49.41 -48.63 41 >50 -28.28 -30.55 -28.64 14 >50 42* *GLIDE docking did not return back any pose for the compound Compound

ACCEPTED MANUSCRIPT Table S5. MM-GBSA binding free energies (kcal mol-1) and in vitro data for Sirt2 GB-8 md -50.04 -54.41 -57.07 -51.15 -58.49 -54.50 -51.78 -51.66 -51.74 -51.80 -55.48 -41.48 -49.74 -54.00 -44.42 -50.09 -46.69 -47.47 -41.42 -51.56 -46.70 -44.98 -47.01 -47.12 -49.33 -44.25 -37.03 -45.12 -47.52 -33.74 -54.34 -33.14 -31.11

md-final -46.74 -56.72 -59.23 -51.10 -62.13 -55.39 -51.50 -51.32 -51.70 -51.84 -56.60 -40.83 -49.56 -54.38 -45.05 -50.00 -48.62 -49.28 -41.09 -51.69 -46.08 -44.66 -47.01 -47.01 -46.90 -46.29 -39.52 -45.81 -48.14 -35.26 -56.93 -32.70 -29.62

RI PT

min -53.64 -57.05 -59.22 -53.25 -60.83 -56.55 -54.09 -53.91 -53.45 -53.98 -57.86 -43.77 -51.28 -56.14 -46.92 -52.50 -49.51 -50.28 -42.32 -53.55 -48.74 -46.62 -48.68 -49.90 -51.86 -47.06 -38.22 -47.24 -49.59 -34.82 -57.93 -33.81 -33.08

SC

md-final -64.51 -67.37 -69.06 -62.56 -68.73 -70.19 -65.25 -62.35 -62.94 -68.83 -69.14 -54.86 -60.77 -63.67 -57.69 -59.25 -58.72 -59.16 -52.79 -64.73 -60.06 -56.49 -56.41 -58.06 -57.79 -57.15 -49.01 -56.78 -57.90 -43.16 -61.15 -43.90 -36.38

M AN U

min -49.15 -59.16 -52.70 -55.98 -61.60 -57.96 -57.18 -55.31 -43.09 -57.11 -42.59 -42.64 -51.63 -47.55 -32.83 -53.18 -48.35 -37.23 -43.47 -40.84 -51.01 -39.91 -35.47 -46.46 -39.12 -35.73 -41.00 -44.20 -29.14 -34.01 -45.62 -29.72 -28.81

GB-1 md -62.88 -69.07 -69.47 -63.19 -68.34 -69.43 -64.71 -62.67 -62.97 -68.40 -66.68 -53.89 -60.38 -62.17 -57.30 -59.89 -60.11 -58.81 -53.01 -64.15 -59.73 -56.25 -59.13 -59.02 -59.79 -55.57 -47.30 -56.20 -59.57 -43.59 -61.43 -43.60 -37.96

TE D

AC C

31 11a 17 11c 18 11d 15a 19 28 11b 29 30 20 37 35 25 32 34 16a 39 26 36 24 33 38 27 21 23 41 22 40 42 14

Activity IC50 (µM) pIC50 8.96 0.0011 8.89 0.0013 8.74 0.0018 8.57 0.0027 8.36 0.0044 8.32 0.0048 8.27 0.0054 8.19 0.0065 8.00 0.0100 8.00 0.0100 7.80 0.0160 7.77 0.0170 7.64 0.0230 7.55 0.0280 7.54 0.0290 7.52 0.0300 7.28 0.0530 7.17 0.0670 6.96 0.1100 6.89 0.1300 6.89 0.1300 6.82 0.1500 6.77 0.1700 6.66 0.2200 6.52 0.3000 6.13 0.7400 5.85 1.4000 5.55 2.8000 5.00 10.0000 4.31 49.0000 >50 >50 >50 -

EP

Compound

ACCEPTED MANUSCRIPT Table S6. Summary of the statistics obtained for the MM-GBSA model based on the lowest binding free energy scores. Sirt3 0.47 (n=28) 0.66 MM-GBSA 31 = 1.99 23 = 2.09 36 = 2.11 38 = 2.15 R2 0.62 (n=24) MM-GBSA* RMSE 0.48 2 q LOO 0.56 * Statistical values calculated for the MM-GBSA method after omitting compounds stated as

RI PT

Statistics R2 RMSE Outliers (XZ-score)

SC

Method

AC C

EP

TE D

M AN U

outliers based on their XZ-score. n; number of compounds used to analyse the statistics.

M AN U

SC

RI PT

ACCEPTED MANUSCRIPT

Figure S6. Principal component analysis (PCA) plot obtained for the studied sirtuin inhibitors by using simple 1D descriptors. The first three principal components, PCA1, PCA2 and PCA3, were projected on the X-, Y- and Z-axis, respectively. Five different clusters of molecules were observed on the

AC C

EP

PCA plot retrieved upon a turn by 45o to the right, c).

TE D

PCA plot, named as G1, G2, G3, G4 and G5 and circled with black. The PCA plot retrieved upon a turn by 45o to the left, a). The PCA plot, b). The