Abstract: Deoxyguanosine kinase is the enzyme catalyzing the first rate limiting step in the salvage pathway of purine nucleotide biosynthesis. Deoxy guanosine 5’ triphosphate, (dGTP) is the axial product of the enzyme reaction. dGTP at increased concentration acts as strong inhibitor of the enzyme. dGTP also inhibits the de novo biosynthesis of purine nucleotides and consequently diminishes the cellular reservoirs of purine nucleotides. Nucleotide insufficiency induced by dGTP prevents DNA synthesis and cellular proliferation. In this way dGTP could be considered for treatment of viral infectious or even cancerous disorders. Molecular Dynamic (MD) study of deoxyguanosine kinase interaction with its inhibitors helps to survey of the structural determinant elements to design more effective inhibitors. In the present work the well known inhibitors of deoxyguanosine kinase including: deoxyguanosine (mono di and triphosphate) and deoxyinosine (mono and tri phosphate) nucleotides for a molecular dynamic simulation was used. Enzyme-Inhibitor complex were obtained using HEX docking software. The best resulting complexes from docking experiments were used for MD simulations. Enzyme-inhibitors complexes were simulated in aqueous solutions at 37°C and 1 atmosphere of pressure. The out put trajectories were used to extract the simulation data. Present results show that nucleotide inhibitors based on their binding pattern could be divided into two classes, deoxyguanosine mono and di phosphates, bind in the proximity of enzyme active site acting as competitive inhibitors and deoxyinosine mono and triphosphate and guanosine triphosphate bind in sites far from the enzyme active site and acting as non competitive inhibitors. Considering present findings it can conclude first, that dGTP the more potent inhibitor is attached to a non-competitive binding site at the C-terminus of enzyme. dGTP is interacting with the negative and positive groups on its binding site, pulling the positive and pushing negative groups of that site deforming the enzyme structure and inhibits the enzyme more effectively. Second, the more negative charges on dGTP in contrast to other inhibitors seems to promote the inhibition more pronounced.
INTRODUCTION
Mammalian deoxyguanosine kinase (dGK), EC 2.7.1.113, with systematic name of ATP: deoxyguanosine 5'-phosphotransferase, is catalyzing the following irreversible reaction:
ATP+deoxyguanosine→ADP+dGMP |
This reaction is the first rate limiting reaction in salvage pathway of purine deoxy nucleotides biosynthesis from the relevant nucleosides. dGK with optimum, pH of 7.5 and 37°C temperature of uses deoxyguanosine as the primary substrate and ATP as phosphate donor for catalysis. dGK is show broad substrate specificity. It uses deoxy adenosine (dA), deoxy cytidine (dC) and deoxy inosine (dI) nucleosides as alternative substrate. dGK is also capable of phosphorylating certain synthetic nucleosides such as ara-G (9-β-D-arabinofuranosyl guanine) and 2-chloro-2'-deoxyadenosine and giving them anticancer or antiviral properties (Kunos et al., 2011; Hanchard et al., 2011; Johansson et al., 2001; Balzarini et al., 2003; Saada, 2008). dGK is constitutively expressed in cellular mitochondrion throughout the cell cycle (Eriksson et al., 2002; Jullig and Eriksson, 2000; Taanman et al., 2003). In human, unlike other organisms, dGK activity is seen predominantly in thymus tissue. Thymus is a specialized gland deals with the processing and maturating T-cells of the immune system. There are miscellaneous reports showing that inherited deficiency of dGK is linked to mitochondrial DNA depletion syndrome. These reports have showing that about 20-30% decrease in dGK expression is sufficient to induce mitochondrial DNA depletion syndrome (Saada-Reisch, 2004; Mandel et al., 2001; Saada et al., 2001; Franco et al., 2007). In catabolic pathway of purine nucleotides, Purine Nucleoside Phosphorylase (PNP), catalyzes the hydrolysis of purine nucleosides to ribose-1-phosphate and free purine bases. In PNP deficient individuals, the cease in nucleotide hydrolysis results in accumulation of deoxy nucleosides e.g., dG. Deoxyguanosine in turn is converted to deoxyguanosine mono phosphate (dGMP) by dGK activity. The resulting dGMP is then converted to dGDP and ultimately to dGTP by successive phosphorylation reactions. The axial product of dGK, i.e., dGTP is a potent inhibitor for the enzyme called ribonucleotide reductase (EC 1.17.4.1). Ribonucleotide reductase is the enzyme responsible for reduction of ribonucleotide to deoxynucleotide in de novo pathway of deoxy nucleotide biosynthesis. Inhibition of de novo pathway by dGTP, consequently unloads the cellular concentrations of deoxynucleotides. Cells with high rate of proliferation such as T cells of the immune system are more susceptible to nucleotide depletion in contrast slowly dividing cells. Accordingly high concentrations of dGTP leads to a dysfunction in T-cell mediated immune response (Lambe et al., 1995; Osborne and Scott, 1983; Leanza et al., 2010). The literature also convey that dGTP is a strong inhibitor for dGK, the first enzyme of salvage pathway of purine biosynthesis, dGTP in addition to the de novo pathway, blocks the salvage pathway as well. Other purine nucleotide including, dGMP, dGDP, dIMP and dITP also show inhibitory effect on dGK activity. These inhibitors are expected to exert the same effect on nucleotide biosynthesis. From medical point of view, there seems to be a paradox regarding the consequences of altered dGK activity. From one prospect, increased activity of dGK helps in activating chemotherapic drugs used for viral infectious or cancer treatment. From the second prospect, dGK inhibition seems to reduce the cellular reservoir of purine nucleotide necessary for DNA synthesis, restricting cellular proliferation. Based on this background dGK inhibitors may become valuable target for more detailed studies, theoretically and experimentally. MD is a powerful technique which is used to simulate the complexes dGK with its inhibitors. MD technique gives us valuable reliable information about enzyme-inhibitor interactions and helps us in understanding the nature and mechanisms of inhibitor binding to dGK (Lengauer and Rarey, 1996; Kitchen et al., 2004; Timmers et al., 2009). The inhibitor constant or Ki, is the concentration of inhibitor required to decrease the maximal rate of the reaction to half of the uninhibited value. Alternatively Ki is the dissociation constant for the enzyme-inhibitor complex or [EI] according to the following equation:
The Ki quantity could be used as a measure for binding tightness of inhibitor to enzyme. Therefore, the lower the Ki of an inhibitor leads to the lower concentration of inhibitor needed to inhibit enzyme activity. The Ki values for dGK were obtained from The comprehensive enzyme information system of BRENDA (www.brenda-enzymes.org) Ki = 3x10-5 mM for dGTP, Ki = 10x10-5 mM for dITP, Ki = 52x10-5 mM for dGDP, Ki = 580x10-5 mM for dGMP and Ki = 3100x10-5 mM for dIMP inhibitors. The objective of the present paper was to interpret experimental Ki of each inhibitor in contrast to other inhibitors according to its binding pattern, structural alteration in enzyme and its energy of binding through expensive calculations of Molecular Dynamic (MD) (Bagchi and Ghosh, 2007).
MATERIALS AND METHODS
Crystal structure of human dGK, PDB ID: 2OCP, with resolution of 2.8 A° and R-value of 0.265 was obtained from the protein data bank, (http://www.rcsb.org/pdb). Coordinate files for nucleotide inhibitors were constructed using ArgusLab 4.0.1 sotware (Mark A. Thompson, Planaria Software LLC, Seattle, WA, http://www.arguslab.com) and energy minimized prior to docking experiments. Docking of nucleotide inhibitors to dGK were carried out using Hex software version 5.1 (Noorbatcha et al., 2009; Smith and Plazas, 2011; Shakyawar et al., 2011). Hex is an interactive program for flexible docking of protein-protein and small ligands to protein either to fixed binding sites or whole the protein structure (Ritchie and Kemp, 2000; Macindoe et al., 2010; Bagchi and Ghosh, 2006). Docking results were scored according to the binding energy and the best solutions of docking were chosen and used for further studies. Nucleotide inhibitors topology and coordinates in GROMACS format were obtained using The Dundee PRODRG2 Server (Schuttelkopf and van Aalten, 2004). Each simulation was initiated with 1500 steps steepest decent minimization followed by 100 ps of system equilibration. A total 10 nanoseconds and for more confirmation 20 nanoseconds simulation was carried out for each of dGK-inhibitors complexes described above at 310 K and 1 atmosphere of pressure. The time steps of 1 femtosecond were applied to all simulations. MD simulations were performed using double-precision MPI version of GROMACS 3.3.1 installed on UBUNTU version 9.10. Each structure of dGK-inhibitors were placed in the center of a cubic box with dimensions of 5.155x6.067x6.189 nm and filled with SPCE water molecules (Maftouni et al., 2011; Bouarkat et al., 2010). LINCS algorithm was used to apply constrain on bonds lengths. The SETTLE algorithm was used to constrain the geometry of water molecules. In the MD protocol, all hydrogen atoms, ions and water molecules were first subjected to 1500 steps of energy minimization by steepest descent and the energy of system were minimized to at least 300 kJ mol-1. The systems then submitted to a short molecular dynamic with position restrains for period of 100 ps and afterwards performed a full molecular dynamics without restrains (Mohammadpour et al., 2011; Zhuang et al., 2011).
RESULTS AND DISCUSSION
As mentioned in introduction section, deoxy purine nucleotides including dGTP, dITP, dGDP, dGMP and dIMP show inhibitory properties on dGK enzyme. dGTP with lower Ki of Ki = 3x10-5 mM is the strongest and dIMP with higher Ki value of Ki = 3100x10-5 mM is the weaker inhibitor among inhibitors. In the present study optimized and well equilibrated structure of dGK were used as a target for docking the inhibitors. Docking results were scored upon their binding energy. The best enzyme-inhibitor complexes were equilibrated for up to 100 ps MD simulation. The final structures of these complexes were used to extract inhibitor binding sites using ArguLab software. The binding site residues are shown in Table 1.
As it depicted in Table 1 all inhibitors did not bind to the same binding sites. The same binding sites of dGDP and dGMP show high residue similarities (in 8 of the total 10-12 residues). This binding site indeed resembles the active site dGK as reported before (Zhu et al., 1998; Eriksson et al., 2001). This finding indicates that dGDP and dGMP are competitive inhibitors for dGK because they bind the active site of dGK. The binding sites of dIMP and dITP also show residues similarity to each other but since their binding site is far from the active site they could be considered as non competitive inhibitors for dGK. The binding site of dGTP is placed in the vicinity of C-terminus and al amino acid. This finding also indicates that dGTP should act as non competitive inhibitor for dGK. The chemical structures of inhibitors reveals that there are some structural differences between inhibitors that may cause them to behave differently. First of all, triphosphate nucleotides (dGTP and dITP) are carrying four negative charges, diphosphate nucleosides (dGDP) are carrying three and mononucleotides (dIMP and dGMP) are carrying only two negative charges. The more negative charges of inhibitors seems to build up a more repulsive force on enzyme through their interaction with negative charges of the binding site. These repulsive forces destabilize the native structure and leads to enzyme inhibition.
Figure 1 shows the changes in the distance between the nearest negative groups in binding sites of inhibitors and nucleotide inhibitors during dynamic simulation. As it evident the distance is increased during simulation depends on the number of negative charges of inhibitors. Molecular dynamic calculations of forces based on Newton first Law are depending on the forces applied and the masses of system particles, therefore we expect the molecular weight of inhibitors should influence the binding energy of enzyme with its inhibitors. Accordingly as we expect the heavier inhibitor, dGTP, shows the stronger inhibition. Following up the change in kinetic energy/total energy ratio during the simulation is used as a criterion for the prevalence of equilibration state. In our system this ration is constant throughout the simulation (data not shown).
| Table 1: | Binding sites residues for dGK inhibitors extracted from last frame of MD trajectories |
| Fig. 1: | Changes in the distance between dGK inhibitors and nearest negative groups on their binding sites |
| Fig. 2(a): | RMSD changes of inhibitors in contrast to dGK protein during simulation obtained for 37°C |
The change in Root Mean Square Displacement (RMSD) for dGK-inhibitor complex during MD experiment is another valuable tool showing the progress of simulation and attainment the stability in the system (Sham et al., 2001). Figure 2a show the change in RMSD of dGK-Inhibitor complexes. This figure show: first, the dGK-inhibitors complex reaches the equilibrium state over 700 ps of simulation, second the more increase in RMSD for inhibitors means the more interaction of inhibitors with dGK during simulation and this is reflected in more changes in dGK to inhibitors distance. Figure 2b shows the change in inhibitors RMSD in contrast to their original places displacement. As shown in Fig. 2b the more potent inhibitor show more displacement from its initial position which in turn means the more penetration of inhibitor in to dGK protein. However dGMP and dIMP show the lower displacements conveying the restricted movement of inhibitors or their entrapment in potential hole inside dGK protein.
Free energy of binding of inhibitors to dGK is calculated based on the method of Linear Interaction Energy (LIE) as reported before by Asi et al. (2004). This kind of binding free energy includes Lennerd-Jones (LJ) and electrostatic contribution.
| Fig. 2(b): | RMSD changes of inhibitors in contrast to their initial position during simulation at 37°C |
| Fig. 3: | Binding free energy of inhibitors to enzyme and its changes through simulation at 37°C |
As depicted in Fig. 3 the binding free energy is increased with increasing inhibitor potency. This finding reveals that dGTP is bind more effectively to dGK and the dGMP and dIMP inhibitors bind weakly in contrast to other inhibitors. Survey of energy contributions in complex formation between inhibitors and dGK shows that more than 95% of binding is belong to electrostatic energy. This finding led us to study the distance changes between inhibitors and the nearest charged residues (positive or negative residues) that present in the binding site, during simulation.
Figure 4 shows the distance changes between dGTP and the nearest residues in binding site. The distance between the negative charge dGTP and the positive charge of Arg3 is decreased with time showing the increased attractive force between the dGTP and dGK. On the other hand the distance between negative charges of dGTP and negative charges of Glu221 is increased with time showing increased repulsive force between dGTP and dGK during simulation at this position of active site. This finding confirm that dGTP deform the enzyme by applying two opposite forces, an attractive force between dGTP and Arg3 and a repulsive force between dGTP and Glu221. These two opposite force may act in enzyme inhibition.
| Fig. 4: | Changes in the distance between dGTP and negative and positive charges in binding site of enzyme |
| Fig. 5: | Mean square displacement changes during simulation for different inhibitors |
Mean square displacement or MSD of inhibitor is the next parameters used in the present study which may be helpful in this context. MSD itself is the average of the distance traveled by inhibitors during simulation and is defined by the following equation (Bamdad et al., 2006):
ri (t)-ri (0) is the distance traveled by inhibitor over simulation time intervals. If the inhibitor encountered no difficulties and penetrate freely in enzyme protein, MSD will increase more quadratically than the restricted ones.
Figure 5 shows the MSD of inhibitors throughout dGK enzyme. As predicted from RMSD changes of inhibitors in contrast to its initial situation, dGTP shows more quadratically diffusing behavior with greatest MDS increase. The rate of MSD increase for each inhibitor with time, is depends on, how often the molecule collide with other particles, so at high density circumstances or in the case of entrapment of inhibitor in to a potential cavity, inhibitor take longer time to diffuse. As we hypothesize for dGMP and its lowest MSD curve show, dGMP probable entrapped in a potential hole and could not get ride of it. Self-diffusion coefficient for inhibitors could be obtained from MSD and using Einstein equations.
As shown in Table 2, the maximum diffusion coefficient of (0.493±0.055)x10-5 cm2 sec were seen for dGTP with freely diffusing properties and the minimum diffusion coefficient of (0.0631±0.029)x10-5 cm2 sec for dGMP with blockage properties. However, other inhibitors as expected behave between these to borders. The next parameter in trajectory analysis of MD simulation is radial distribution function or RDF for inhibitors. RDF is a statistical mechanics parameter obtained by MD simulation, describes how the atomic density of inhibitors varies as a function of the distance from dGK protein during the simulation. The obtained results are shown in Fig. 6 for inhibitors accordingly.
The lower height or lower distribution of the curve at the maximum point means the more distribution of inhibitor through protein and more penetration of inhibitor into dGK. From this point of view, dGTP show more going through protein behavior than other inhibitors and the rest of inhibitors show degreased penetration behavior as expected definitely. Based on RMSD, MSD and radial distribution changes upon interactions of inhibitors with dGK we expect to see structural alterations in protein structures reflected whether, in secondary and/or tertiary structures. The first parameter used to study the tertiary structure changes upon dGK interaction with nucleotide inhibitors is Accessible Surface Area (ASA) which is the surface area in Å2 of dGK protein that is accessible to a non aqueous solvent (Lee and Richards, 1971). Any change in ASA during simulation indicates structural changes in tertiary structure. Since most ASA curves obtained in MD simulation for proteins are broad and constitute of overlapping bands it is hard to analyze them simply and extract reliable results even for comparative purpose. Using different smoothing methods may helps to pickup better curves for comparison but the process may be misleading. We found that the variance in ASA changes during the simulation gives a better parameter to show the extent of inhibitor interaction with protein. Seemingly the more variance in ASA belongs to the more interacting inhibitor through tertiary structure.
Table 3 show the variance in ASA during simulation for inhibitors. dGTP, the stronger inhibitor pose the greatest variance in ASA and increase ASA to about 2.53%.
| Table 2: | Diffusion coefficient for dGK inhibitors |
| Table 3: | Variance of ASA in contrast to its initial value in percent |
| Fig. 6: | Radial distribution curve for inhibitors across the simulation box at 37°C |
The increase in ASA is in agreement with conformational changes in dGK from native toward slightly unfolded structure with inhibited activity. dGMP and dIMP as could be expected show the lowest variance in ASA which in turn agree with their weak inhibitory effect. In the cases of dGDP and dITP the variance in ASA are not changed as anticipated. It is of great importance to note that the mechanism of inhibition of dGK may differ from one inhibitor to another inhibitor, e.g., some inhibitors may act via tertiary, some via secondary and some via both structures simultaneously, so dGDP and dITP may exert their inhibitory effect on secondary and tertiary structures differently. dGK structure shows that this enzyme in its secondary structure contains, 9 alpha helix comprising about 51.09% of protein sequence, 3 beta strands comprising 12.27% and random coil of about 36.64% of total protein residues.
Table 4 shows the changes of protein secondary structure percentage during simulation. As shown in table 4 dGMP shows no effect on secondary structures of dGK and no changes were seen, while dIMP affects the alpha helix structure, decreased it to about 5.13% in expense of increasing random coil to about 7.14%. Even though dGTP slightly increase random coil (about 1.19%) in expense of decreasing alpha structures but dGTP primarily impose its inhibitory effects via tertiary structure distortion. Hydrogen bond change during MD simulation is another important factor to study in this work.
Table 5 summarizes hydrogen bond variance between protein-protein and protein-solvent along simulation in the presence of inhibitors. As it evident the most variance in protein-protein and protein-solvent hydrogen bonds is seen for dGTP which is affects both secondary and tertiary structures in synergistic fashion. However the pattern of hydrogen bonds change is dependent to how the inhibitors distorted the protein structure. The radius of gyration of a protein or Rg is the root mean square distance of protein atoms from protein centroid. Rg describes the overall spread of the molecule and determines the protein structure compactness, the more Rg the less compactness of protein structure (Lee and Richards, 1971). The Rg change curve of dGK during the MD simulation has also broad and overlapping essence so we obliged to use Rg variance instead of simple Rg for comparative purpose.
Table 6 accumulates the Rg variance for dGK in the presence of inhibitors. dGTP and dITP the two potent inhibitors show the more variance in Rg during the simulation.
| Table 4: | Secondary structure changes in the presence of inhibitors (%) |
| Table 5: | Variance of hydrogen bonds between protein-protein and protein-solvent in the presence of inhibitors |
| Table 6: | Percent of variance in gyration radius of dGK in the presence of inhibitors |
CONCLUSION
From the present findings, it can be concluded taht, dGTp the more potent inhibiter is attached to a non-competitive binding site at the C-terminus of enzyme. The dGTp is interacting with the negative and positive groups on its binding site, pulling the positive and pushing negative groups of that site deforming the enzyme structure and inhibits the enzyme more effectively. Moreover, the more negative charges on dGTp in contrast to other inhibitors seems to promote the inhibition more pronounced.