Physical nature of intermolecular interactions inside Sir2 homolog active site: molecular dynamics and ab initio study

In the present study, we analyze the interactions of NAD+-dependent deacetylase (Sir2 homolog yeast Hst2) with carba-nicotinamide-adenine-dinucleotide (ADP-HPD). For the Sir2 homolog, a yeast Hst2 docking procedure was applied. The structure of the protein–ADP-HPD complex obtained during the docking procedure was used as a starting point for molecular dynamics simulation. The intermolecular interaction energy partitioning was performed for protein–ADP-HPD complex resulting from molecular dynamics simulation. The analysis was performed for ADP-HPD and 15 amino acids forming a deacetylase binding pocket. Although the results indicate that the first-order electrostatic interaction energy is substantial, the presence of multiple hydrogen bonds in investigated complexes can lead to significant value of induction component.


Introduction
Weak noncovalent interactions involving molecules of biological importance have been identified to play an important role in many fundamental processes in nature. The detailed analysis of strength of hydrogen bonds and London dispersion interactions can provide a better understanding of the structure and stability of biomolecules [1]. The forces are responsible for specificity of DNA-protein binding. Many pharmaceutical ligand-protein interactions are noncovalent in their nature. Although the X-ray measurements results are often available and give a lot of geometrical information, the magnitude and character of interaction between the substrate-specific ligand and enzymes need further exploration.
Due to their biochemical properties Sir2 (silent information regulator 2) proteins are involved in various biological processes such as DNA metabolism, regulation, and repair of double-stranded breaks [2,3]. This group of proteins can be involved in transcriptional silencing, apoptosis and chromosome stability. Their promoting activity in longevity in yeast mother cells was also proved [4]. It has been also indicated that some Sir2 homologs in eukaryotes have been implicated in the proper cell cycle progression, radiation resistance, and genomic stability [5]. To date, investigations strongly suggest that their biological activity is largely dependent on their deacetylate properties and is also partially modulated by concentration of nicotinamide in cells. The conformation of Sir2 is unique among all Sir proteins and was carefully determined [6]. In the current study we present the structural and energetical consequences of interactions between NAD+-dependent deacetylase (Sir2) with carba-nicotinamide-adeninedinucleotide (ADP-HPD) [7]. According to Sanders' study, nicotinamide molecule is bound to the D-pocket which has hydrophobic character due to the PHE67, PHE184, around the pyridine ring and PHE44 together with ILE117 proximal to the carboxamide moiety (See Fig. 1). The docking procedure and molecular dynamics simulations were performed in order to gain an insight into the structural and energetical basis of Sir2 enzyme inhibition. Moreover, the total intermolecular energy between protein and ADP-HPD was divided into physically relevant terms determining the nature of interactions in complex in question.

Docking procedure
During the docking procedure, the structure of the Sir2 homolog with an acetyllysine histone H4 peptide yeast Hst2 was used. The structural data of the considered protein, denoted as 1SZC [8], was downloaded from the Protein Data Bank. Both structures used during the docking procedure, namely Sir2 protein and optimized ADP-HPD ligand, contained only polar hydrogen atoms. The docking stage was performed with use of the united-atom scoring function implemented in AutoDock Vina [9]. The space declared during the docking procedure (18 × 20 × 24) has been limited to the active site of Sir2 enzyme. The docking procedure for considered subunits was repeated ten times with exhaustiveness equal 30; its outcomes are presented in Fig. 2. In all cases, atoms of ADP-HPD ligand backbone are localized almost in this same point of conformational space of the active site. Small differences are mainly related to localization of hydrogen atoms.

Molecular dynamics simulation
The structure of the Hst2-ADP-HPD complex obtained during docking procedure was used as a starting point for molecular dynamics simulation. Ligand structure was characterized with the use of amber force field parameters and the atomic charges were calculated according to the Merz-Kollman scheme via the RESP procedure [10] at HF/6-31G* level. The complex structure was neutralized and immersed in a periodic box of a TIP3P water box [11]. The analyzed system was heated to 300 K during 100 ps of MD simulation and the temperature was controlled by a Langevin thermostat [12]. The periodic boundary conditions and SHAKE algorithm [13] were applied for 60 ns of molecular dynamic simulation. The analysis of molecular dynamics simulation parameters, namely potential and kinetic energies of system and RMSD values, confirm that the first 20 ns of obtained trajectory is an equilibration period. Values presented in Fig. 3 confirm that both elements of the considered system, namely protein and ADP-HPD ligand, after this stage reached dynamic stability. The final 40 ns of trajectory, collected during molecular dynamics simulation, was used in analysis of the interaction between considered subunits. Structural analysis was performed with use of the VMD package [14]. In all molecular dynamics simulations, the AMBER 11 package was used [15].

Intermolecular interaction energy partitioning
Intermolecular interaction energy partitioning was performed for the protein-ADP-HPD complex resulting from molecular dynamics simulation. The analysis was performed for ADP-HPD and amino acids forming a deacetylase binding pocket. In the present study, we discuss the results obtained using the variational-perturbational scheme implemented in a modified version of the GAMESS US package [16,17].
In this approach, the total intermolecular interaction energy calculated at Hartree-Fock level of theory (ΔE HF ) is decomposed into the first-order electrostatic (ε el (10) ), Heitler-London exchange (ΔE ex HL ) and delocalization components (ΔE del HF ). The delocalization term, also referred to as the deformation component, contains charge-transfer and polarization effects. Since the dimer-centered basis set is employed during the calculations, both the total intermolecular interaction energy as well as its components are free from the basis set superposition error.
The analysis was performed for ADP-HPD and 15 amino acids forming a deacetylase binding pocket. The RHF/aug-cc-pVDZ level of theory and about 1700 contracted Cartesian Gaussian basis functions were employed to describe the systems. Due to the size of the studied complexes, we have not attempted to determine interaction energy components resulting from the partition of the correlation energy at the MP2 level of theory (i.e., uncoupled dispersion energy, correlation correction to electrostatic energy and exchange-delocalization term) [18].   medium and weak interactions. The next non-negligible hydrogen bonds are created by SER270 and THR37 but in both cases frequency of occurring and binding power is much smaller than in the case of ASN248 and LEU249. The middle orthophosphoric part of ADP-HPD is stabilized by interactions, of strong and middle binding power, involving ALA33. Presented distribution of hydrogen bonds indicate on small mobility of ADP-HPD in conformational space of active site what is confirmed also by stable RMSD values (av. 1.23 ± 0.18). The next group of important forces related to stabilization of considered ligand is hydrophobic interactions. Active site of Sir 2 contains many amino acids, which could be involved in such interactions. The most important role can be fulfilled by aromatic amino acids localized in the vicinity of heterocyclic rings of ADP-HPD, namely phenyloalanies PHE44, PHE67, and PHE184. Figure 5 presents distances between aromatic systems of amino acids and nicotinamide part of considered ligand during molecular dynamics simulations. Presented distributions indicate that at least one of considered phenylalanine moieties plays an important role in stabilization of the ADP-HPD molecule. In the case of PHE44 by 50 % of simulation time, the distance between aromatic rings does not exceed 5 Å, and mutual orientations between considered molecules allow the occurrence of stacking interactions. Also, localization of PHE67 in the active site during simulation ensures a large group of conformations, enabling the occurrence of hydrophobic interactions with heterocyclic ring of ADP-HPD.

Selection criteria
During this stage we selected the structure of the complex of ADP-HPD with Sir2 protein, which should be the most representative relative to the conformations occurring during considered simulation times. For this purpose, we used structural data describing distributions of 11 of the most important hydrogen bonds and distances between aromatic rings of the considered subunits. Several criteria were considered to choose a representative structure. One of the conditions that should fulfill the chosen structure was the smallest total value of considered hydrogen bonds, which should potentially correlate with the highest affinity of ligand toward active site. The second criteria required that hydrogen bond length and also distance between the considered aromatic subunits should belong to the most-represented intervals of distributions obtained during simulation times.

Intermolecular interaction energy partitioning
The existence of multiple hydrogen bonds between the ADP-HPD molecule and the surrounding environment suggests that electrostatics should be important. Indeed, in terms of absolute numbers, its value is quite large, and has mainly a stabilizing character. Its highest value was found for the positively charged ARG45-ADP-HDP complex. As might be expected, the lowest values of electrostatic component (close to 0 kcal/ mol) were observed for complexes formed by the aromatic The data are given in kcal/mol ring of phenylalanine and ADP-HDP molecule. As is shown in Table 2, in the case of some complexes, the stabilization due to the mutual polarization of interacting subsystems reflected in the negative value of delocalization component (encompassing induction terms) and is also quite important.
In the case of complexes containing ALA33 and ASP118, stabilization originates strongly from the interaction of permanent and induced multipole moments of subsystems due to their mutual polarization. In the case of phenylalanine and arginine complexes, the effects resulting from quantummechanical exchange of electrons between subsystems are particularly small. Some results indicate that although the first-order electrostatic is substantial, it might be completely quenched by associated Pauli repulsion. It was observed, for example, for weak interactions of ASN245 and TYR269 with ADP-HDP molecules. Consequently, it would be interesting to estimate the correlation contributions to the total intermolecular stabilization.
In order to characterize the relationships between intermolecular interaction energy components, Spearman rank correlation was used [19].
The data collected in Table 3 were obtained using the expression above, where d is the difference between ranks of the compared terms and n is the number of investigated complexes. It is seen that among of all pairs of components, the correlation coefficients of ΔEel/ΔEdel and ΔEel/ΔEHF are high. Between the exchange and delocalization components, an evident negative correlation with coefficient equal to −0.8 is observed (Table 4).
The nature of the hydrogen bonds stabilizing of investigated system can be characterized based on the theory of intermolecular interactions. Since the higher-order delocalization effects involve charge transfer excitation, they can be associated with covalent character of interactions. It was also shown by Grabowski et al. that delocalization has a dominant stabilizing character in the group of complexes with short dihydrogen bonds (distances close to 1 Å) [20]. As they observed for the remaining hydrogen bonds, electrostatic energy should be the most important. Based on the mentioned study, we can conclude that the hydrogen bonds that are partially covalent have an Edel/Eel ratio higher than 0.45. In the case of our systems, consequently, the most covalent characters have ALA33 and SER225 complexes.

Conclusions
Each fragment of ADP-HPD ligand, namely nicotinamide, adenine nucleotides, and orthophosphoric chains, play important roles in stabilization of a complex with enzyme by numerous interactions of varying nature. The presented values describing time evolution of ADP-HPD molecule in the complex with Sir2 protein correlate with outcomes of docking procedure. Most all of the important interactions identified on the docking stage are present during molecular dynamics simulations. The greatest accumulation of interactions is observed in the space of D-pocket due to the presence of numerous hydrogen bonds (ILE117, ASP118, ALY16 B) and hydrophobic interactions (PHE44, PHE67), also in other regions of active site are important bindings present by most simulation time. Quite small and stable RMSD values of ADP-HPD molecule, supported by high frequency of bindings creation during whole simulation, indicate small mobility and limited conformational flexibility of considered ligands in the space of enzyme active site.
As is seen, even at the uncorrelated level of theory, the picture of interactions is quite complicated. Although in the case of the studied complex typical stacking motifs are not present between ADP-HPD molecules and neighboring amino acids, it follows from previous studies that one might expect non-negligible additional stabilization resulting from dispersion component of interaction energy [21]. The analysis of correlations between the intermolecular interaction components revealed existing dependence between positive and negative for ΔEel/ΔEdel and ΔEex/ΔEdel, respectively.