|Home | About | Journals | Submit | Contact Us | Français|
Following our previous study (Piquemal et al., New J. Chem., 2003, 27, 909), we present here a DFT study of the inhibition of the Tyrosinase enzyme. Broken-symmetry DFT computations are supplemented with Spin-Flip TD-DFT calculations, which, for the first time, are applied to such a dicopper enzyme. The chosen biomimetic model encompasses a dioxygen molecule, two Cu(II) cations, and six imidazole rings. The docking energy of a natural substrate, namely phenolate, together with those of several inhibitor and non-inhibitor compounds, are reported and show the ability of the model to rank the most potent inhibitors in agreement with experimental data. With respect to broken-symmetry calculations, the Spin-Flip TD-DFT approach reinforces the possibility for theory to point out potent inhibitors: the need for the deprotonation of the substrates, natural or inhibitors, is now clearly established. Moreover, Electron Localization Function (ELF) topological analysis computations are used to deeply track the particular electronic distribution of the Cu-O-Cu three-center bonds involved in the enzymatic Cu2O2 metallic core (Piquemal and Pilmé, J. Mol. Struct.: Theochem, 2006, 77, 764). It is shown that such bonds exhibit very resilient out-of-plane density expansions that play a key role in docking interactions: their 3D-orientation could be the topological electronic signature of oxygen activation within such systems.
Tyrosinase is a monophenol oxygenase (IUBMB enzyme Nomenclature: EC 18.104.22.168) found widely within the living realm (Lerch, 1987; Dugas, 1989; Solomon et al., 2001). This enzyme is involved, for instance, in the formation of the melanin shield protecting plants from external stress, and is also responsible for the brown color of fruits upon contact with air (O2). Human Tyrosinase can be found in skin melanocytes or in the iris of eyes; it is responsible for the dark color of the cerebellum. Several mutations are associated with known pathologies, such as albinism which is linked to a thermal instability of the enzyme at physiological temperature. [Giebel et al., 1991 and references therein]. The chemical role of Tyrosinase is to catalyze the ortho-hydroxylation of monophenols into ortho-diphenols or, alternatively, the oxidation of ortho-diphenols into the corresponding ortho-quinones. In mammals, it performs for instance the hydroxylation of L-tyrosine into 3-(3,4-dihydroxyphenyl)-L-alanine (DOPA) which is oxidized into DOPA quinone (catechol oxidase activity), as depicted in Figure 1 (Decker and Tuczek, 2000; Kim and Uyama, 2005; and references therein).
From a structural viewpoint, tyrosinase is classified as a type III copper enzyme because of the presence of two coupled copper cations in its active site, whose role is to activate dioxygen to initiate catalytic activity.(Ross et al., 1991). Two forms of the enzyme have thus to be considered: a deoxy state for the native protein, and an oxy state when dioxygen binds. The deoxy form is EPR-silent, and its UV-visible spectra does not exhibit any charge transfer bands. These findings are consistent with a singlet spin state and formal Cu(I)/3d10 redox states of the copper cations. The oxy form also is EPR-silent, but exhibits a much more complex electronic structure: it encompasses a peroxide ligand (O22−) bridging two formal Cu(II)/(3d9) cations exhibiting a strong antiferromagnetic coupling (Gherman and Cramer, 2009; Piquemal et al., 2006; and references therein).
The structure-function analysis of the enzyme has been hampered for several years due to the lack of any crystallographic data: the first crystallographic structure has been reported only recently (Matoba et al., 2006). Nevertheless, several cross data including genetic sequence homology and X-ray spectroscopy (X absorption, EXAFS and XANES) had shown strong structural similarities between Tyrosinase and Hemocyanin (Decker and Tuczek, 2000 and references herein). The latter is also a class III dicopper protein and plays a role of oxygen carrier in mollusk and arthropod hemolymph. When the structure of Panulirus Interruptus hemocyamin was solved (Volbeda and Hol, 1989), it allowed a clear-cut picture of the bimetallic active site that appeared to involve six histidine residues. Based on these structures several bio-inspired models have been synthesized (Kitajima and Morooka, 1994; Karlin et al., 1998; Murthy et al., 2001; Palavinici et al., 2005; Mirica et al., 2006; Tolman, 2006; and references therein), some of them have been found to hydroxylate phenol derivatives to the corresponding quinones. Extensive theoretical studies were also performed focusing on the structure and reactivity of compounds exhibiting a Cu2O2 core. The computational approaches however have to face the complex electronic structures associated with such bimetallic cores and adequate quantum chemistry tools need to be employed (Gherman and Cramer, 2009). However, beyond the interest of theoreticians for such subtle electronic structures and as the crystal structure has been recently solved (Decker et al., 2006), little remains known about the effective mechanism of the enzyme. Indeed, understanding and inhibiting Tyrosinase would be important in medicine due to its clear involvement in Parkinson’s disease (Xu et al., 1997), melanoma (Prezioso et al., 1992), and hyperpigmentation phenomenon (Maeda and Fukuda, 1991). Moreover, inhibiting Tyrosinase can prevent the unwanted darkening of fruits and seafood, which has an important financial cost (Friedman, 1996). A computational approach would thus help to understand the details of such inhibition
Building on the immense amount of work available, we started to work on this problem in 2003 (Piquemal et al., 2003) and showed on the basis of Density Functional Theory (DFT) calculations, that it is possible to predict computationally the inhibition of a model of the enzyme by 2-aminophenol (2-APOH). It was first shown that both the substrate and the inhibitor should be deprotonated to form a stable complex with the active site. Second, it was shown that only phenolate binds to the oxy form. We finally suggested a competitive inhibition mechanism relying on the deprotonation of the substrates. Recently, such predictions have found a clear experimental confirmation (Mirica et al., 2006). In the present contribution, we report an extension of our previous study to other inhibitors, and compare the theoretical inhibition hierarchy to the experimental one. To improve our previous work, two levels of computations are employed here: a broken-symmetry approach (BS-DFT) and the recently introduced spin-flip TD-DFT technique (Shao et al., 2003; Wang and Ziegler, 2004; Neese, 2009), hereafter referred to as SF-DFT. To go deeper into the understanding of the inhibition mechanism, the analysis of the electronic structures is complemented by the Electron Localization Function (ELF) topological analysis (Becke and Edgecombe, 1990; Silvi and Savin, 1994). Such an analysis is expected to provide some insight into the particular role of the peculiar 3-center Cu-O-Cu bonds within the Cu2O2 core (Piquemal and Pilmé, 2006) in the regioselective approach of the substrate or inhibitors.
Several inhibitors having demonstrated experimental potency thanks to their structural similarity to phenol, a model of the tyrosine side-chain, will be tested under their deprotonated form, namely 2-aminophenol (2-APOH), benzoic and o-toluic acids; kojic acid (5-hydroxy-2-hydroxymethyl-4-pyrone) used in cosmetics, and 4-hydroxybenzaldehyde. A known non-inhibitor compound, 2,6-dimethylbenzoic acid, is also investigated to test the ability of the model to discriminate between inhibitor and non-inhibitor molecules. All compounds are depicted on Figure 2.
The paper is organized as follows: we will first present the model of the tyrosinase active site, describe the used methodologies, and then discuss the results.
The bio-inspired model of the oxy-tyrosinase active site, hereafter referred to as TYR-O2, encompasses two copper cations, six imidazole ligands, and dioxygen: the resulting complex [(Cu(ImH)3)2-O2]2+ recovers that of our previous study (Piquemal et al., 2003). Such a structure was initially retained by Lind et al. (1999) to model the six histidine residues of the active site. From a structural viewpoint, this model has been shown to be in fine agreement with experimental data and is consistent with the recent X-ray structure of tyrosinase (Decker et al. 2006). As previous studies showed that the (μ-O)2 arrangement is not catalytically relevant (Tolman, 2006), we only consider here the μ2:η2-η2 form. The substrate and the inhibitors, taken in their deprotonated forms, are presented on Figure 2. No solvent effects were considered in the previous investigations as the active site of tyrosinase is known to be encapsulated in a hydrophobic pocket of the whole protein (Conrad et al., 1994).
Complete geometry optimizations were carried out using Jaguar (Jaguar 4.1, 2000) at the DFT level of computation using the B3LYP functional (Becke, 1988; Lee et al., 1998) with the LACVP** basis set (McLean et al., 1980; Rasolov et al., 1998). Broken symmetry solutions were sought to account for the open-shell nature of the Cu2O2 singlet core. Indeed, despite its known spin delocalization problem, BS-DFT has been shown to be able to capture the main electronic features of type III copper active sites (Piquemal and Pilmé, 2006; Gherman et al., 2009; and references therein). In our BS-DFT computations, spin density values on the oxygen and copper atoms are consistent with those reported for comparable active site models (Cramer et al., 2006). Actually, recent studies even showed that DFT performs better than expected to reproduce accurate structures (Rode and Werner, 2005; Cramer et al., 2006) when compared to high-level multireferential post-Hartree-Fock approaches. It has moreover been shown by means of the ELF topological analysis that the essential electronic features of the central Cu2O2 core, especially the three-centre nature of the Cu-O-Cu bonds, were conserved when going from a CASSCF to a DFT formalism (Piquemal and Pilmé, 2006).
Zero-point vibrational energy (ZPE) and basis set superposition error (BSSE) are included in the reported energies. ZPE and entropies were evaluated within the harmonic approximation at T = 298 K.
To go a step further towards accuracy we have also considered the recent Spin-Flip Time-Dependent DFT formalism (SF-DFT) as an alternative to the BS-DFT techniques (de la Lande et al., 2007; de la Lande et al., 2008a; de la Lande et al., 2008b; de la Lande et al., 2009). This approach was shown to bevery promising for the treatment of mononuclear Cu/O2 adducts found in non-coupled copper monooxygenases. The present study thus offers the opportunity to test this formalism in the case of binuclear sites. All SF-DFT computations were performed at the BS-DFT minima using the Q-Chem package (Shao et al., 2006).
In order to follow the organization of the electron density at the Cu2O2 core by means of ELF topological analyses, we have used a modified TopMod package (Noury et al., 1999; Pilmé and Piquemal, 2008). We here briefly recall that the ELF function can be interpreted as the signature of the electronic-pair distribution. Its topological analysis allows the partition of space into domains, the so-called “basins”, which are related to the intuitive chemical picture of bonds. For a given atom A, core electrons generate core basins labelled C(A). Regions associated with lone pairs are labelled V(A) and bonding regions denoting chemical bonds between atoms A and B are noted V(A,B), or V(A,B,C) in the case of three-center bonds. These ELF regions match closely the domains of the VSEPR model. Details about ELF and its applications to bioinorganic chemistry can be found in a recent review paper (Piquemal et al., 2008).
We computed the interaction energies (ΔE) between the tyrosinase model and its substrate/inhibitor:
To take thermodynamics (and therefore for entropic effects) into account, ΔE energies include several corrections: ΔΔG computations are then provided (at 298K). Final ΔΔG values include the ZPE and BSSE energies added to the inhibitor deprotonation energy denoted as Edeprot as the substrate and the inhibitor should be deprotonated to form a stable complex with the active site (Piquemal et al., 2003). Energetic results are gathered in Tables 1 and and22.
2-APO− clearly appears as the best inhibitor. Its stabilization is 20 kcal/mol larger compared to that of the phenolate natural substrate. 2-APO− is followed by deprotonated kojic acid and 4-hydroxybenzaldehyde, the latter being favoured by its weak deprotonation energy. It is important to point out that SF-DFT discriminates more clearly 2-APO− from the other inhibitors: an additional −12.4 kcal/mol stabilization is found whereas deprotonated kojic acid is disfavoured by +9.4 kcal/mol.
Benzoate and o-toluate appear as the worst inhibitors. Toluate presents an energetic balance slightly superior to phenolate by 0.2 kcal/mol at the BS-DFT level and by −0.2 kcal/mol with SF-DFT. The inclusion of thermal corrections switches their relative stability: deprotonated o-toluic acid is now stabilized by −0.5 kcal/mol. However, considering the levels of calculations and their intrinsinc uncertainty coupled to thermal agitation under physiological conditions, we can consider that these two compounds have the same affinity for the active site.
Dimethylbenzoate does not appear to inhibit the active site model: it has a highly unfavorable energetic balance compared to phenolate.
In all the structures, the Cu2O2 cores exhibit a strong butterfly effect oriented towards the inhibitors as can be seen on Figures 3 and and4.4. Detailed structural data are given in Table 3 for the optimized geometries. The more pronounced effect is observed for the most potent inhibitor, i.e. 2APO−: the Cu-Cu bond length (4.31Å) is increased when compared to that in the phenolate complex (3.33Å). It can be noted that the oxygen atoms of the deprotonated kojic acid appear in interaction with the hydrogen atoms of the imidazole residues.
The optimized structure of the dimethylbenzoate complex shows a Cu2O2 metallic core presenting a deformation which is visually closer to a “zig-zag” form than to the usual μ–η2:η2 binding mode: the carboxylate oxygen atoms are oriented towards the copper cations. It is important to note that, despite the elongated O-O distance, an ELF analysis shows a conservation of the topological organization of the electronic structure of the μ–η2:η2 binding mode described previously (Piquemal and Pilmé. 2006). This topology involves two three-center Cu-O-Cu bonds being above and below the Cu2O2 central metallic core. This striking feature of the Cu2O2 core displayed in Figure 6 appears to be conserved in all the inhibitor/active site model calculations. In all cases, the Cu-O-Cu bond basin is oriented towards the anionic oxygen atoms of the inhibitor pointing out the importance of this unusual local organization of the electronic densities. For the non-inhibitor compound, (dimethylbenzoate complex), carboxylate oxygen atoms are oriented towards the copper atom and not toward the Cu-O-Cu bond. Nevertheless, the topology involving two three-center Cu-O-Cu bonds above and below the Cu2O2 central metallic core is conserved despite the apparent “zig-zag” configuration.
Such a “resilient” topological organization of the electron density (“out-of-plane” vs.” in-plane”) of dioxygen is likely to be linked to the “oxygen activation” and, then, may directly influence the regioselectivity of the enzyme by governing the orientation of the substrate/inhibitor approach.
Overall, the results are satisfactory and the model is able to qualitatively discriminate between inhibitor and non inhibitor compounds (Table 4). We now address the quantitative agreement of the computational results with the available experimental data (Maddaluno and Faull, 1988; Conrad et al., 1994; Bubacco et al., 2000; Jimenez et al., 2001).
Several apparent inhibition constants (Ki) have been reported for Tyrosinase in the literature (Table 5). A note of caution has however to be brought as the genes encoding for Tyrosinase sometimes come from different living species. Here we consider the measurements reported for mushroom Tyrosinase (Maddaluno and Faull, 1988; Conrad et al., 1994; Jimenez et al., 2001) and for Streptomyces antibioticus Tyrosinase (Bubacco et al. 2000). This is a critical point: even if the active site is highly conserved between different species, the inhibition activity of a given drug may differ, probably due to the immediate proteic environment. For example, results for the benzoic acid indicate a good inhibiting power for the mushroom but a lesser one for Streptomyces Antibioticus. This observation is in agreement with a study which associated its inhibition capability to its ability to interact favourably with the hydrophobic enzymatic pocket (Conrad et al., 1994).
The role of entropy has also to be taken into account. However, our model encompasses only a limited number of atoms through an harmonic approximation. Indeed, as can be seen from Table 6 (experimental results for the mushroom), o-toluate appears to be a more efficient inhibitor compared to benzoate thanks to a more negative ΔS value and despite a less favourable ΔH value, which finally results in an overall more favorable free enthalpy (ΔG = ΔH − T. ΔS). This trend is also observed for Streptomyces Antibioticus: o-toluate is by far a better inhibitor than benzoate. From these results, and considering that our calculations do not take into account the stabilization (or destabilization) contribution that could be brought to the energy by the environment, we have to expect an inversion of benzoate and o-toluate.
The inhibition ranking obtained from our calculations (Table 4) then appears in relatively good agreement with experiment. First, the two inhibitors known to be really potent, namely 2-APO− and deprotonated kojic acid, are ranked in the two first places, which is satisfactory. The increased theoretical stabilization of 2-APO− obtained at the SF-DFT level compared to the BS-DFT approach is more in line with the extremely small experimental Ki constant. Second, the non-inhibitor compound dimethylbenzoate appears last in the ranking: this shows the capability of the model to discriminate between an inhibitor and a compound without action on the enzyme.
Finally, deprotonated 4-hydroxybenzaldehyde appears to be slightly more stabilized than benzoate probably because of the absence of lateral chains in our model. Indeed, these latter have been shown to stabilize the benzoate adduct (Conrad et al., 1994).
Our study of the tyrosinase shows that it is possible to use quantum chemistry to perform in silico predictive docking computation using a relatively small active site model compared to reality. SF-DFT computations were for the first time performed on such an enzyme model. Results on enzyme/substrate interactions confirm that BS-DFT is able to reasonably describe such systems as a first theoretical inhibition hierarchy has been computed in reasonable agreement with experimental inhibition constants. SF-DFT computations were shown to be able to more clearly discriminate the best inhibitor.
The present computations tend to confirm that reactivity occurs when the active site is in oxygenated form and that the substrate and inhibitors are active in their deprotonated form, which confirms our previous hypothesis. This study also illustrates the importance of taking into account environment and entropy effects, and more generally the need to consider a greater number of residues to refine the model if aiming for a finer virtual screening of the inhibitors.
The ELF topological analysis has shown that the metallic core conserves a non conventional three-center bond organization along the docking process even in the case where a zig-zag configuration seems to be present. This very resilient topology clearly indicates that this peculiar density organization plays a direct role in the interaction as the Cu2O2 core acts as a single entity. As dioxygen is known to be activated in these complexes, such Cu2O2 three-center bonds could be the topological signature of dioxygen activation. Another aspect that has to be investigated is the efficiency of alternative spin states in the oxidative reaction of the substrate. This work is currently in progress and will be reported in due time. Thanks to the availability of a crystal structure, future work will focus on the inclusion of the whole protein in the simulations using QM/MM techniques and free energy perturbation to include environment dynamics and entropy. Such an approach should enable the screening of dozens of compounds, and maybe even more as the available computational resources increase.
These results will also be important for the design of advanced molecular modelling techniques able to handle the difficult electronic structure of the enzyme (Piquemal, J.-P., Williams-Hubbard, 2003; Gresh et al., 2007; Diedrich et al. 2008).
The computations were performed at CRIHAN (Saint-Etienne-du-Rouvray, France) on project 2008011, IDRIS (Orsay, France) and NIEHS (RTP, NC, USA). This research was supported by the Intramural Research program of the NIH and NIEHS (USA). The authors are deeply indebted to Dr Claude Giessner-Prettre who initiated the project.