Comparative Assessment of Theoretical Methods for the

Jun 29, 2007 - Globally, B3LYP had the best average performance in the test, closely followed by MP2, whereas B1B95 was the least accurate method. The...
5 downloads 12 Views 149KB Size
9146

J. Phys. Chem. B 2007, 111, 9146-9152

Comparative Assessment of Theoretical Methods for the Determination of Geometrical Properties in Biological Zinc Complexes Se´ rgio Filipe Sousa, Pedro Alexandrino Fernandes, and Maria Joa˜ o Ramos* REQUIMTE, Departamento de Quı´mica, Faculdade de Cieˆ ncias, UniVersidade do Porto, Rua do Campo Alegre, 687, 4169-007 Porto, Portugal ReceiVed: April 1, 2007; In Final Form: May 16, 2007

In the present study, we have compared the performance of the density functional theory (DFT) functionals B1B95, B3LYP, B97-2, BP86, and BPW91 with MP2 for geometry determination in biological mononuclear Zn complexes. A total of 15 different basis sets, of rather diverse complexity, were tested, several which included also three different types of common effective-core potentials: Los Alamos, Steven-Basch-Krauss, and Stuttgart-Dresden. In addition, the ability to describe mononuclear Zn biological systems using relatively simple models of the metal coordination sphere, comprising only the metal atom and a simplified representation of the ligands at the first coordination sphere, starting from a set of high-resolution X-ray crystallographic structures, is evaluated for 90 combinations of method/basis set. The results show that the use of such models allows for a relatively accurate description of the Zn-ligand bond lengths, although failing to correctly represent the topology of the metal coordination sphere (namely, the angles involving the metal atom) if constraints at the CR atoms are not considered. Globally, B3LYP had the best average performance in the test, closely followed by MP2, whereas B1B95 was the least accurate method. The study also points out B3LYP/CEP121G and B3LYP/SDD, which use, respectively, the Steven-Basch-Krauss and the Stuttgart-Dresden effective-core potentials, as the best compromise between accuracy and CPU time for the geometrical characterization of metal-ligand bond lengths in Zn biological systems.

Introduction Zinc is a biologically essential transition element and an integral component of more than 400 enzymes responsible for a multitude of functions, present in virtually all forms of life.1-3 In fact, zinc, the second most abundant element in biology, is the only metal known to have representatives in each of the six fundamental classes of enzymes established by the International Union of Biochemistry.2 Approximately 2800 human proteins have been identified to bind zinc in vivo by a recent bioinformatics survey, representing something like 10% of the total human proteome.4 Carboxypeptidase A, thermolysin, carbonic anhydrase, and alcohol dehydrogenase are just some of the most well-studied Zn enzymes, but the list of zinc metalloenzymes and metalloproteins is long and extremely diversified.3,5 Zinc has some very specific properties that differentiate it from the other first-row transition metals and that make zinc a very attractive target for biological systems, explaining its almost ubiquitous presence in organisms. Among these is the fast exchange of ligands, the lack of redox activity, the flexible coordination geometry, the role as a Lewis acid, and the notably high bioavailability.3,6 In fact, under physiological conditions, this metal exists only in the dicationic state, which has a closed d shell (d10 configuration). The Zn(II) ion is thus diamagnetic and forms colorless complexes. Furthermore, no naturally occurring Zn(II) complexes with (colored) tetrapyrrole ligands exist. These characteristics rather limit the range of experimental techniques directly applicable to the study of Zn biological systems and make Zn one of the “spectroscopically quiet” metals in biology.7 It is therefore not surprising that the application of * To whom correspondence should be addressed. E-mail: mjramos@ fc.up.pt.

theoretical and computational methods to the study of Zn biological systems has flourished so early in the history of the field.8-11 Over time, several pungent questions regarding specific biological Zn systems, both at the structural and mechanistic levels, were unraveled only by application of computational methodologies. Recent examples include alcohol dehydrogenase,12 matrix metaloproteinase,13 Cu,Zn-superoxide dismutase,14,15 peptide deformylase,16 and the puzzling enzyme farnesyltransferase.17-20 When studying reaction pathways, the search for the intermediates and transition states requires a large number of calculations to be performed. Hence, it is usual to keep the size of the model used as small as possible to allow a comprehensive analysis of the aspects under study21 taking into consideration that the choice of the model must include all essential groups to ensure that the chemistry of the system is accurately represented. This requirement makes the study of biological systems that involve metal atoms particularly difficult, as the treatment of the metal (or metals, in the case of polynuclear centers) significantly increases the computational expense associated. In metalloproteins, the electronic structure of the metal atom, together with its immediate environment, is normally responsible for the lion’s share of the catalytic process.21 A detailed representation of the metal and of its first coordination sphere is therefore, in many cases, sufficient to adequately represent the reactions taking place at the enzyme’s active site. The inclusion of specific second-shell residues, whose influence in the reaction is known from experiment or is suspected, is also a standard procedure in the computational treatment of metal sites. Until now, most of the comprehensive studies involving biological Zn systems have relied almost essentially on this type

10.1021/jp072538y CCC: $37.00 © 2007 American Chemical Society Published on Web 06/29/2007

Theoretical Methods for Biological Zn Centers of models, the so-called cluster models. An alternative approach to the use of cluster models in the treatment of Zn biological systems is the use of hybrid quantum mechanics/molecular mechanics (QM/MM) methods, which are able to ally to the accuracy of a QM method, the low computational cost characteristic of the molecular mechanical methodologies. However, even though QM/MM calculations can nowadays be routinely applied to the study of enzymes comprised simply by standard amino acid residues, the application of QM/MM calculations to the study of biological systems involving one or more transition-metal atoms is still, in many cases, extremely demanding and cumbersome.21 Parameters for the metal atom and the directly coordinating residues need to be derived and added to the molecular mechanical force field. These parameters typically involve the metal bond lengths and angles and the correspondent force constants and partial charges and are normally specific for each enzyme, or even for each enzymatic state, particularly when the identity of the ligands at the metal coordination sphere changes during catalysis.22-25 Despite the existing difficulties, some studies using QM/MM-based methodologies in transition-metal biochemistry have been reported, albeit with some limitations, in the literature.26,27 The use of small active site models treated with density functional theory (DFT)-based methodologies is still paramount in the computational study of biological metal systems21,28 being routinely applied in the study of enzyme active sites, in the determination and validation of enzymatic reaction mechanisms, in the resolution of mechanistic paradoxes and dilemmas, in the characterization of metal-dependent structural and spectroscopic properties, and even in the determination molecular parameters for the application of molecular mechanical methods. In the particular case of Zn biological systems, the d10 configuration of the metal atom results in a ligand field stabilization energy of zero for all possible geometries, leading to spatially not directed (isotropic) polarization effects. Therefore, no geometry is intrinsically more stable than another, which confers to zinc metalloenzymes an ability to accommodate several different coordination environments during an enzymatic catalytic process without significant energetic cost. Biological Zn complexes are hence typically characterized by a high flexibility of the metal coordination sphere, the lack of a specially preferred geometrical arrangement by the metal atom, particularly flat potential energy surfaces, and a multitude of low-lying energy states connected by low-energy barriers. These aspects render the determination of the correct geometry, for both minima and transition states, a particularly thorny task. Moreover, the resulting structures are typically highly dependent on the combination of method/basis set used. The performance of several DFT functionals has been recently reviewed and compared in detail in the literature, in terms of binding energies, geometries, and barriers heights. Most of these studies focus exclusively on main group elements.29-35 Some efforts in transition-metal chemistry have also been reported, but these studies have relied almost essentially in metal dimers.36-39 This study compares the application of five different DFT functionals (B1B95, B3LYP, B97-2, BP86, and BPW91) and MP2 theory in combination with 15 different basis sets ranging from a minimal basis set to more computational demanding triple-ζ basis sets and passing through three different types of common effective-core potentials (Los Alamos,40-42 StevenBasch-Krauss,43-45 and Stuttgart-Dresden46,47). Rather than performing a benchmarking analysis of the performance of the several methods and basis sets considered, this study further aims to evaluate the ability to describe mononuclear Zn

J. Phys. Chem. B, Vol. 111, No. 30, 2007 9147 TABLE 1: Summary of the DFT Functionals Testeda functionals

year

type

χ

exchange functional

correlation functional

B1B95 B3LYP B97-2 BP86 BPW91

1996 1994 2001 1988 1991

HM-DFT H-DFT H-DFT pure pure

25 20 21 0 0

Becke88 Becke88 B97-2 Becke88 Becke88

Becke95 Lee-Yang-Parr B97-2 Perdew86 Perdew-Wang91

a HM, hybrid meta; H, hybrid; χ, percentage of HF exchange in the functional.

biological systems using relatively simple models of the metal coordination sphere, comprising only the metal atom and a simplified representation of the ligands at the first coordination sphere, starting from a selection of high-resolution X-ray crystallographic structures. In particular, the average errors introduced by the use of such models with the 90 combinations of method/basis set in the determination of Zn-ligand bond lengths and angles are compared and discussed. Finally, comments on the validity of such an approach for the characterization of the metal coordination sphere on Zn biological systems are drawn. Computational Methods Calculations were carried out using the Gaussian 03 suite of programs.48 Five different commonly available DFT functionals were tested (Table 1). The performance of MP2 calculations in the treatment of the type of systems considered in this study was also evaluated. The DFT functionals tested were B1B95,49,50 B3LYP,49,51-53 B97-2,54,57 BP86,49,55 and BPW91.49,56 B1B95 is a one-parameter hybrid meta-DFT functional that incorporates the B88 exchange functional of Becke49 and the B95 correlation functional.50 B3LYP is by far the most popular hybrid DFT functional and was designed as a three-parameter hybrid functional51 that employs the B88 exchange functional of Becke49 and the Lee, Yang, and Parr correlation functional.52 B97-2 uses Becke’s 97 hybrid exchange-correlation energy functional57 as it was parametrized by Wilson et al.54 BP86 is a pure DFT functional, which uses Becke’s 1988 gradient-corrected exchange functional (B)49 together with Perdew’s 1986 gradient-corrected correlation functional (P86).55 Finally, BPW91 is a pure DFT functional that combines Becke’s 1988 gradient-corrected exchange functional (B)49 with Perdew-Wang’s 1991 correlation functional.56 A total of 15 different basis sets were assessed in this study. These are listed in Table 2 together with the valence basis set considered for each of the atom types that make up the Zn complexes in the set of biological models studied. This number includes six basis sets that use effective-core potentials (ECPs). These are the CEP-4G,43-45 CEP-31G,43-45 CEP-121G,43-45 LanL2MB,40-42 LanL2DZ,40-42,58 and SDD.46,47,58 The basic principle behind ECP calculations is that the chemically inert core electrons, whose description typically requires a large set of gaussians, can be substituted by an approximate function, the so-called pseudopotential or effective-core potential. This process significantly reduces the number of electrons that are explicitly described in the calculation, thereby lowering the computational expense associated. In fact, many chemical properties such as bond strengths, electron affinities, polarizabilities, ionization potentials, and molecular geometries are essentially determined by the valence electrons. The use of ECPs has thus established itself as a particularly useful tool in the treatment of systems involving heavy atoms and most notably in the case of transition metals and is generally considered a

9148 J. Phys. Chem. B, Vol. 111, No. 30, 2007

Sousa et al.

TABLE 2: Description of the Several Basis Set Considered valence electrons (if ECPs used) basis set

basis set type

STO-3G 3-21G 6-31G 6-31G(d) 6-311G(d,p) CEP-4G CEP-31G CEP-121G LanL2MB LanL2DZ SDD SV SVP TZV TZVP

all electrons all electrons all electrons all electrons all electrons small core small core small core large core large core small core all electrons all electrons all electrons all electrons

Zn

20 20 20 12 12 20

S

6 6 6 6 6

C, N, O

4, 5, 6 4, 5, 6 4, 5, 6

valence basis set Zn

S

C, N, O

H

3/33/3/3 3/33/21/21 6/6631/31 6/6631/31/1 6/6631/31/1 4211/411 4211/411 4211/411 3/2/5 21/11/41 311111/22111/411 63311/53/41 63311/53/41/1 842111/631/411 842111/631/411/1

3/33 3/321 6/631 6/631/1 631111/42111/1 4 31 121 3/3 21/21 531111/4211 5311/511 5311/511/1 73211/6111 73211/6111/1

3/3 3/21 6/31 6/31/1 6/6311/1 4 31 121 3/3 721/41 6111/41 511/31 511/31/1 62111/411 62111/411/1

3 21 31 31 311/1 3 31 311 3 31 31 31 31/1 311 311/1

TABLE 3: Summary of the Biological Zn Systems Tested PDB Code

class

enzyme

resolution (Å)

Zn coordination sphere

refs

1HET 2TLX 1M4L 2CBA 1D8W 1EVL 1VHH

I. oxidoreductase I. oxidoreductase III. hydrolase IV. lyase V. isomerase VI. ligase signaling protein

alcohol dehydrogenase thermolysin carboxypeptidase A carbonic anhydrase II L-rhamnose isomerase threonyl-tRNA synthethase sonic hedgehog

1.15 1.65 1.25 1.54 1.60 1.55 1.70

Cys Cys-His-Wa His His-Glu-W His His-Glu-W His His-His-W Glu-Asp-His-As p His His-Cys-W His His-Asp-W

60 61 62 63 64 65 66

a

W, a water molecule.

good and safe choice when properties related to the valence electron system are to be investigated. The basis sets tested in this study include three different types of effective-core potentials. CEP-4G, CEP-31G, and CEP-121G use the StevenBasch-Krauss pseudopotentials43-45 (sometimes credited as SKBJ), LanL2MB and LanL2DZ use the Los Alamos pseudopotentials40-42 (sometimes referred as Hay and Wadt), whereas SDD uses the Stuttgart-Dresden pseudopotentials46,47 (also known as Stoll-Preuss or simply SP). Details concerning these pseudopotentials and the corresponding schemes included in Gaussian 03 for the treatment of Zn and of the directly coordinating atoms are described in Table 2 together with a comparison with the several all-electron basis sets considered in this study. No diffuse functions were considered in any of the basis sets of this study. Although in the calculation of reaction energies, barrier heights, electron affinities, ionization potentials, and related properties the use of diffuse functions can normally yield better results, for the determination of bond lengths and angles the use of such functions leads typically to a rather marginal improvement on the results with a significant increase on the computational cost.59 Accordingly, in this study, the attention was focused on the effect of the size, type, and contraction scheme of the basis set and on the inclusion of polarization functions. Models for the mononuclear biological Zn systems were prepared from a selection of high-resolution crystallographic structures chosen from the ones with the best resolution within the Protein Data Bank (PDB). Care was taken to ensure that this selection included examples from several different classes of Zn enzymes and that an appropriately diverse set of coordination spheres was considered. Details on the biological Zn systems chosen are presented in Table 3. Figure 1 illustrates the models considered. Starting from the initial PDB structures, models including only Zn and the first coordination sphere of the metal atom were prepared. Conventional modeling of the amino acid side chains was used, that is, the Zn ligands aspartate

(or glutamate), cysteinate, and histidine were modeled by acetate, methylthiolate, and imidazole, respectively. Each of these initial structures was then energy minimized using B1B95, B3LYP, B97-2, BP86, BPW91, and MP2 together with each of the 15 basis sets presented in Table 2. The resulting structures were compared with the initial models, and the metal-ligand bond lengths and angles were compared with the values reported in the correspondent high-resolution X-ray structures. Results and Discussion Table 4 presents the mean signed errors (MSEs) calculated for each of the 90 combinations of method/basis set, considering all the Zn-ligand bond lengths in the seven models of Zn biological systems considered in this study, which were prepared from seven high-resolution X-ray crystallographic structures selected from the PDB. The MSE is taken as the difference between the values calculated with the combination method/ basis set for the model systems and the correspondent values reported in the initial X-ray structures. A negative MSE indicates that the application of a given combination method/basis set to the type of models considered underestimates the value (in this case, the bond length), whereas a positive MSE indicates that the value is overestimated. Table 5 exhibits the correspondent mean unsigned errors (MUE). Globally, the results show that the use of relatively simplified models of the metal coordination sphere comprising only the metal atom and the directly coordinating ligands to describe Zn biological systems results in fairly accurate metal-ligand bond lengths. In fact, in the set of 6 methods and 15 basis sets tested, the average MUE is 0.16 Å, and the use of a minimal basis set (STO-3G) yields an average MUE of 0.33 Å. We would like to recall that for an X-ray structure obtained at a medium level of resolution (around 2.0 Å), the expected uncertainty in metal bond lengths is at least 0.3 Å.67,68 So, even the use of a minimal basis set can, in principle, yield more accurate metalligand bond lengths than a low resolution (>2.5 or 3 Å) X-ray

Theoretical Methods for Biological Zn Centers

J. Phys. Chem. B, Vol. 111, No. 30, 2007 9149

Figure 1. Models evaluated in this study to describe mononuclear Zn biological systems and correspondent PDB codes.

TABLE 4: Mean Signed Error (MSE) in Zn-Ligand Bond Lengths (Å)

TABLE 5: Mean Unsigned Error (MUE) in Zn-Ligand Bond Lengths (Å)

B1B95 B3LYP B97-2 BP86 BPW91 MP2 average STO-3G 3-21G 6-31G 6-31G(d) 6-311G(d,p) CEP-4G CEP-31G CEP-121G LanL2MB LanL2DZ SDD SV SVP TZV TZVP

-0.34 -0.24 -0.11 -0.17 -0.15 -0.17 -0.09 -0.10 0.01 -0.14 -0.18 -0.10 -0.13 -0.10 -0.26

-0.29 -0.25 -0.07 -0.17 -0.07 -0.05 -0.09 -0.06 -0.09 -0.11 -0.08 -0.09 -0.14 -0.07 -0.05

-0.32 -0.22 -0.09 -0.12 -0.05 -0.16 -0.19 -0.16 0.00 -0.03 -0.16 -0.10 -0.15 -0.07 -0.04

average

-0.14 -0.10 -0.11 -0.10 -0.11 -0.09 -0.12

-0.24 -0.18 -0.15 -0.07 -0.15 -0.08 -0.08 0.03 -0.10 -0.16 -0.08 -0.06 -0.08 -0.05

-0.24 -0.16 -0.14 -0.14 -0.14 -0.14 -0.06 0.03 -0.02 -0.14 -0.09 -0.10 -0.11 -0.13

-0.18 -0.12 -0.05 -0.14 -0.18 -0.17 -0.07 -0.16 -0.07 -0.04 -0.07 -0.06 -0.10 -0.06 0.00

-0.28 -0.22 -0.11 -0.15 -0.11 -0.14 -0.11 -0.10 -0.02 -0.07 -0.13 -0.09 -0.11 -0.08 -0.09

structure. The use of better basis sets can compete in terms of accuracy with medium resolution X-ray structures. Hence, it is not surprising that the application of QM calculations to improve the accuracy of protein X-ray structures has developed as a particularly interesting strategy with huge potential.69 One of the conclusions that can de drawn from the analysis of Table 4 is that the use of the type of models considered in

B1B95 B3LYP B97-2 BP86 BPW91 MP2 average STO-3G 3-21G 6-31G 6-31G(d) 6-311G(d,p) CEP-4G CEP-31G CEP-121G LanL2MB LanL2DZ SDD SV SVP TZV TZVP average

0.37 0.25 0.13 0.21 0.19 0.21 0.10 0.13 0.10 0.17 0.19 0.13 0.18 0.12 0.29 0.17

0.33 0.26 0.10 0.18 0.09 0.14 0.11 0.10 0.18 0.16 0.10 0.12 0.15 0.11 0.09 0.13

0.37 0.24 0.11 0.18 0.11 0.22 0.20 0.17 0.10 0.09 0.19 0.13 0.17 0.11 0.08 0.15

0.26 0.19 0.18 0.09 0.21 0.10 0.14 0.13 0.15 0.17 0.11 0.09 0.10 0.09 0.15

0.25 0.18 0.18 0.17 0.21 0.17 0.10 0.11 0.09 0.17 0.12 0.13 0.13 0.16 0.15

0.24 0.15 0.09 0.19 0.22 0.22 0.09 0.18 0.19 0.08 0.09 0.10 0.17 0.08 0.08 0.14

0.33 0.23 0.13 0.19 0.14 0.20 0.13 0.14 0.14 0.12 0.15 0.12 0.15 0.11 0.13 0.16

this study typically results in underestimated metal-ligand bond lengths. This effect was observed for virtually all combinations of methods and basis sets and is a result of the influence of the enzymatic environment, which forces a distortion of the metal complexes from the ideal geometry, weakening the correspondent bonds, and resulting in an increase on the bond lengths. This effect is well described in the literature for the case of the Zn enzyme farnesyltransferase.25

9150 J. Phys. Chem. B, Vol. 111, No. 30, 2007 One of the problems of considering the MSE in the analysis of the performance of the several combinations of methods and basis sets is that it is frequent in some lower quality approaches with large positive and negative errors for the errors to cancel out yielding a low MSE, whereas higher quality approaches with smaller errors but all of the same sign exhibit higher MSEs. These features are partially observable in Table 4. Hence, despite the usefulness of MSEs to evaluate the bias in the results that arise from the use of small active site models to represent Zn biological systems and whether values are under- or overestimated, the use of MSUs (as presented in Table 5) allows a more coherent analysis of the overall performance of the several different combinations tested. The results presented in Table 5 show that the MUE for the Zn-ligand bond lengths varies from a maximum of 0.37 Å to a minimum of 0.08 Å. B3LYP had the best performance with an average MUE of 0.13 Å in the set of basis sets considered (STO-3G not included) followed by MP2 (0.14 Å) and by BP86, B97-2, and BPW91 (0.15 Å). B1B95 had the worst performance (MUE ) 0.17 Å). BPW91 and B1B95 rendered major optimization problems with some basis sets (particularly for STO-3G and the ones considering ECPs) for some of the complexes (especially 1EVL and 2CBA) and in some cases were unable to find an energy minimum. MP2 was the most computationally demanding method, whereas B3LYP was, by far, the fastest. In terms of the basis sets, some interesting conclusions can also be inferred from Table 5. An analysis of the average MUE for the several basis sets, considering the six methods used, reveals a significant decrease on the average MUE when moving from an STO-3G to a 3-21G basis set and from this basis set to a 6-31G basis set. However, the use of more computationally demanding basis sets like 6-31G(d) or 6-311G(d,p) does not necessarily improve the accuracy on the metal-ligand bond lengths. In fact, the use of 6-31G(d) with B1B95, B3LYP, B972, and MP2 resulted in an important increase in the MUEs over the 6-31G basis set, whereas the use of the 6-311G(d,p) basis set resulted in an improvement over the 6-31G basis set only in the case of B3LYP, BP86, and BPW91. The performance of the Ahlrichs70,71 basis sets (SV, SVP, TZV, and TZVP) was also typically better than 6-31G and 6-311G(d,p). Interestingly, the inclusion of polarization functions (as in SVP and TZVP) did not improve the accuracy in the determination of metalligand bond lengths of the more simple descriptions (SV and TZV). To analyze in more detail these unexpected results, we have compared the performance of all the combinations methods/basis sets for a set of four zinc dimers (Zn-Zn, ZnNe, Zn-Ar, and Zn-Kr)72-75 and with the results obtained by Zhao and Truhlar37 using 18 density functionals with the augcc-PVTZ basis set for the same species. The same general effect was observed in these small molecules when moving for more complete basis sets. Details are presented in the Supporting Information section. A particularly interesting result, for the field of metalloproteins, was the performance of the several ECP basis sets in the study of the seven models of Zn biological complexes. In fact, the ECP basis sets used were able to obtain results of at least the same level of accuracy of the best all-electron basis sets considered in this study. With the exception of the CEP-4G, the remaining basis sets considering ECPs resulted in an average MUE between 0.12 and 0.15 Å. The LanL2DZ basis set, which in the formulation included for Zn in the Gaussian 03 program is a large core basis set and which includes pseudopotentials also in the elements from the second period (sulfur), achieved the best result in this regard (average MUE of 0.12 Å).

Sousa et al. TABLE 6: Ranking of the 90 Combinations of Method/Basis Set Tested on the Basis of the Mean Unsigned Error (MUE) Calculated for the Zn-Ligand Bond Lengths (Å)a B1B95 B3LYP B97-2 BP86 BPW91 MP2 average STO-3G 3-21G 6-31G 6-31G(d) 6-311G(d,p) CEP-4G CEP-31G CEP-121G LanL2MB LanL2DZ SDD SV SVP TZV TZVP average a

88 82 37 73 65 75 24 35 23 55 69 38 62 33 85

86 84 18 60 14 42 29 16 64 47 17 32 45 25 9

87 79 28 58 30 77 71 53 22 5 66 36 56 27 1

89 83 70 59 11 74 21 43 39 46 54 31 7 19 10

89 81 63 57 52 72 50 20 26 6 49 34 40 41 48

80 44 12 68 76 78 8 61 67 3 13 15 51 2 4

6

1

4

3

5

2

15 14 6 12 9 13 4 8 7 3 11 2 10 1 5

Special emphasis is dedicated to the top 20 combinations.

TABLE 7: Average Unsigned Error in Zn-Ligand Angles (°) B1B95 B3LYP B97-2 BP86 BPW91 MP2 average STO-3G 3-21G 6-31G 6-31G(d) 6-311G(d,p) CEP-4G CEP-31G CEP-121G LanL2MB LanL2DZ SDD SV SVP TZV TZVP average

45.0 46.1 48.6 43.0 41.6 42.5 47.5 47.0 41.1 46.4 45.7 43.4 42.2 48.2 43.6 44.8

42.9 41.1 45.5 43.9 49.4 42.6 48.7 42.1 42.4 46.7 47.4 47.0 46.8 47.8 46.3 45.5

44.7 46.3 47.5 41.9 43.6 42.3 48.8 47.6 42.9 45.8 43.7 48.9 46.1 45.8 45.9 45.5

46.3 47.8 44.2 44.0 42.7 47.4 38.4 41.9 47.0 46.9 49.3 46.9 49.4 46.4 45.6

46.3 46.5 44.1 44.4 42.1 46.3 41.9 42.7 46.2 46.4 49.7 47.8 48.5 46.9 45.7

41.3 44.9 45.8 42.6 43.7 40.5 45.0 42.1 40.7 45.2 43.5 46.4 42.5 48.5 42.3 43.8

43.5 45.2 47.0 43.3 44.5 42.1 47.3 43.2 41.9 46.2 45.6 47.4 45.4 48.0 45.2 45.1

Table 6 presents a ranking of the 90 combinations of method/ basis set tested on the basis of the MUE calculated for the Znligand bond lengths (Å) on the seven Zn biological complexes tested. Although such hierarchy should be taken with care, the results show that MP2 has seven combination method/basis sets in the top 20 of the total 90 combinations tested, whereas B3LYP has five. The best combination involving the B1B95 functional occupies a modest 23rd position (LanL2MB). Taking into consideration the computational expense associated to the several basis sets, and that B3LYP is by far the fastest of the six methods tested, our results suggest B3LYP/CEP-121G and B3LYP/SDD as the best compromise between accuracy and CPU time for the geometrical characterization of metal-ligand bond lengths in Zn biological systems. The overall trends described in detail for the full set of Zn-ligand bond lengths were individually observed also with minor exceptions for the Zn-O, Zn-N, and Zn-S bond lengths. Details are presented in the Supporting Information. Table 7 presents the MUEs calculated for the angles where Zn is the central atom. The results show that in Zn biological systems the use of simple models of the metal coordination sphere, comprising only the metal atom and a simplified representation of the ligands, results in an average MUE on the angles of 45.1°. Moreover, the results show that these high MUEs are virtually independent of the quality of the combination of method/basis set used, with MUE varying between 38.4° and 49.7°. In fact, the lowest MUE was obtained with the

Theoretical Methods for Biological Zn Centers combination BP86/CEP-121G, whereas B3LYP/6-311G(d,p), BP86/TZV, and BPW91/SV rendered the highest MUEs. This evidence highlights one of the limitations of the use of small active site models in the treatment of Zn biological systems. In fact, the structural effect of the enzyme on the topology of the active site Zn coordination sphere is significantly more important (average MUE on the angles of ca. 47%) than its effect on the metal-ligand bond lengths (average MUE less than 6%). Studies of Zn biological systems involving aspects that are highly dependent on the orientation of the ligands at the Zn coordination sphere should take this into consideration. However, this problem can be normally overcome by the use of constraints at the CR atoms (or terminal atoms) of the amino acid residues that are directly coordinated to the metal atom. Conclusions In this paper, we have compared the performance of the DFT functionals B1B95, B3LYP, B97-2, BP86, and BPW91 with MP2 for geometry determination in biological mononuclear Zn complexes using a total of 15 different basis sets. In addition, we have evaluated the ability to structurally describe mononuclear Zn biological systems using relatively simple models of the metal coordination sphere, comprising only the metal atom and a simplified representation of the ligands at the first coordination sphere. The results show that in the set of 90 combinations of methods and basis sets tested, the use of such simplified models allows for relatively accurate description of the metal-ligand bond lengths (average MUE of 0.16 Å) but fails to correctly represent the topology of the metal coordination sphere (average MUE in angles where the metal is the central atom of 45.1°) if additional precautions (for example, constraining the position of the CR) are not taken into consideration. Furthermore, the results show that the use of these models typically results in underestimated metal-ligand bond lengths, in comparison with the real biological system (average -0.12 Å). Globally, B3LYP had the best performance of the test (average MUE of 0.13 Å), closely followed by MP2 (average MUE of 0.14 Å), whereas B1B95 had the worst set of results (average MUE of 0.17 Å), showing that B3LYP, despite its known drawbacks, remains a useful tool in the theoretical characterization of Zn complexes at least in terms of the topology of the system. For other properties, however, particularly in the case of ionization potentials, electron affinities, heats of formation, and most notably barrier heights, the use of more recent and sophisticated density functionals should be taken into account. In particular, BB1K, MPW1K, and MPWB1K have been shown to reproduce barrier heights with an MUE of less than 1.5 kcal/mol (against an MUE of 4.5 kcal/mol by B3LYP) in a data set of 76 reactions of molecules composed of main group elements.35 The good performance of these more recent density functionals is also well portrayed for a variety of properties in several other studies, although transition metals are normally absent from the test systems.32,34,38,39,59,76,77 In these studies, however, the gain in accuracy in the determination of bond lengths and bond angles by using these more sophisticated density functionals is typically much smaller.38,39,59,76-78 In terms of the basis sets tested, the results show that the use of more computationally demanding basis sets like 6-311G(d,p) and TZVP does not necessarily improve the accuracy on the metal-ligand bond lengths over the values typically obtained with a more computationally accessible basis set like 6-31G. Even in the cases where the use of a more computational expensive basis set leads to an increase in the quality of the

J. Phys. Chem. B, Vol. 111, No. 30, 2007 9151 results, the improvement on the metal-ligand bond lengths is normally remarkably small. Furthermore, the use of effectivecore potentials is a particularly interesting approach, which allows a significant gain in terms of CPU time without noticeable loss of accuracy in comparison with the more traditional all-electron basis sets. Taking into consideration the computational expense involved in the application of methods/ basis sets tested, our results suggest B3LYP/CEP-121G and B3LYP/SDD (MUE of 0.10 Å) as the best compromise between accuracy and CPU time for the geometrical characterization of metal-ligand bond lengths in Zn biological systems. Acknowledgment. The authors thank the FCT (Fundac¸ a˜o para a Cieˆncia e a Tecnologia) for financial support (POCI/ QUI/61563/2004). Se´rgio Filipe Sousa also acknowledges the FCT for a doctoral scholarship (SFRH/BD/12848/2003). Supporting Information Available: Results for the ZnZn, Zn-Ne, Zn-Ar, and Zn-Kr dimers. Tables with the individual average unsigned errors in Zn-O, Zn-N, and Zn-S bond lengths (Å) for the models of the Zn biological complexes. This material is available free of charge via the Internet at http:// pubs.acs.org. References and Notes (1) Vallee, B. L.; Auld, D. S. Biochemistry 1990, 29, 5647. (2) Vallee, B. L.; Auld, D. S. Proc. Natl. Acad. Sci. U.S.A. 1990, 87, 220. (3) Lipscomb, W. N.; Strater, N. Chem. ReV. 1996, 96, 2375. (4) Andreini, C.; Banci, L.; Bertini, I.; Rosato, A. J. Proteome Res. 2006, 5, 196. (5) Coleman, J. E. Curr. Opin. Chem. Biol. 1998, 2, 222. (6) McCall, K. A.; Huang, C. C.; Fierke, C. A. J. Nutr. 2000, 130, 1437. (7) Penner-Hahn, J. E. Coord. Chem. ReV. 2005, 249, 161. (8) Bertini, I.; Luchinat, C.; Rosi, M.; Sgamellotti, A.; Tarantelli, F. Inorg. Chem. 1990, 29, 1460. (9) Nakagawa, S.; Umeyama, H. J. Theor. Biol. 1982, 96, 473. (10) Demoulin, D.; Pullman, A. Theor. Chim. Acta 1978, 49, 161. (11) Liang, J. Y.; Lipscomb, W. N. Biochemistry 1987, 26, 5293. (12) Agarwal, P. K.; Webb, S. P.; Hammes-Schiffer, S. J. Am. Chem. Soc. 2000, 122, 4803. (13) Diaz, N.; Suarez, D.; Sordo, T. L. J. Phys. Chem. B 2006, 110, 24222. (14) Pelmenschikov, V.; Siegbahn, P. E. M. Inorg. Chem. 2005, 44, 3311. (15) Amadei, A.; D’Alessandro, M.; Paci, M.; Di Nola, A.; Aschi, M. J. Phys. Chem. B 2006, 110, 7538. (16) Leopoldini, M.; Russo, N.; Toscano, M. J. Phys. Chem. B 2006, 110, 1063. (17) Sousa, S. F.; Fernandes, P. A.; Ramos, M. J. Biophys. J. 2005, 88, 483. (18) Sousa, S. F.; Fernandes, P. A.; Ramos, M. J. J. Comput. Chem. 2007, 28, 1160. (19) Sousa, S. F.; Fernandes, P. A.; Ramos, M. J. J. Mol. Struct. (THEOCHEM) 2005, 729, 125. (20) Sousa, S. F.; Fernandes, P. A.; Ramos, M. J. Proteins 2007, 66, 205. (21) Himo, F. Theor. Chem. Acc. 2006, 116, 232. (22) Comba, P.; Remenyi, R. Coord. Chem. ReV. 2003, 9, 238. (23) Ryde, U. Proteins 1995, 21, 40. (24) Banci, L. Curr. Opin. Chem. Biol. 2003, 7, 143. (25) Sousa, S. F.; Fernandes, P. A.; Ramos, M. J. Theor. Chem. Acc. 2007, 117, 171. (26) Noodleman, L.; Lovell, T.; Han, W. G.; Li, J.; Himo, F. Chem. ReV. 2004, 104, 459. (27) Cavalli, A.; Carloni, P.; Recanatini, M. Chem. ReV. 2006, 106, 3497. (28) Sousa, S. F.; Fernandes, P. A.; Ramos, M. J. J. Am. Chem. Soc. 2007, 129, 1378. (29) Zhao, Y.; Pu, J.; Benjamin, J. L.; Truhlar, D. G. Phys. Chem. Chem. Phys. 2004, 6, 673. (30) Zhao, Y.; Truhlar, D. G. J. Phys. Chem. A 2006, 110, 10478. (31) Zhao, Y.; Truhlar, D. G. J. Chem. Theory Comput. 2006, 2, 1009. (32) Zhao, Y.; Truhlar, D. G. J. Chem. Theory Comput. 2005, 1, 415. (33) Zhao, Y.; Truhlar, D. G. J. Phys. Chem. A 2005, 109, 6624. (34) Zhao, Y.; Truhlar, D. G. J. Phys. Chem. A 2005, 109, 5656.

9152 J. Phys. Chem. B, Vol. 111, No. 30, 2007 (35) Zhao, Y.; Gonzalez-Garcia, N.; Truhlar, D. G. J. Phys. Chem. A 2005, 109, 2012. (36) Zhao, Y.; Truhlar, D. G. J. Phys. Chem. A 2006, 124, 224105. (37) Zhao, Y.; Truhlar, D. G. J. Phys. Chem. A 2006, 110, 5121. (38) Schultz, N. E.; Zhao, Y.; Truhlar, D. G. J. Phys. Chem. A 2005, 109, 11127. (39) Schultz, N. E.; Zhao, Y.; Truhlar, D. G. J. Phys. Chem. A 2005, 109, 4388. (40) Hay, P. J.; Wadt, W. R. J. Chem. Phys. 1985, 82, 299. (41) Wadt, W. R.; Hay, P. J. J. Chem. Phys. 1985, 82, 284. (42) Hay, P. J.; Wadt, W. R. J. Chem. Phys. 1985, 82, 270. (43) Cundari, T. R.; Stevens, W. J. J. Chem. Phys. 1993, 98, 5555. (44) Stevens, W. J.; Krauss, M.; Basch, H.; Jasien, P. G. Can. J. Chem. 1992, 70, 612. (45) Stevens, W. J.; Basch, H.; Krauss, M. J. Chem. Phys. 1984, 81, 6026. (46) Andrae, D.; Haussermann, U.; Dolg, M.; Stoll, H.; Preuss, H. Theor. Chim. Acta 1990, 77, 123. (47) Dolg, M.; Wedig, U.; Stoll, H.; Preuss, H. J. Chem. Phys. 1987, 86, 866. (48) Frisch, M. J.; Trucks, G. W.; Schlegel, H. B.; Scuseria, G. E.; Robb, M. A.; Cheeseman, J. R.; Montgomery, J. A.; Vreven, T.; Kudin, K. N.; Burant, J. C.; Millam, J. M.; Iyengar, S. S.; Tomasi, J.; Barone, V.; Mennucci, B.; Cossi, M.; Scalmani, G.; Rega, N.; Petersson, G. A.; Nakatsuji, H.; Hada, M.; Ehara, M.; Toyota, K.; Fukuda, R.; Hasegawa, J.; Ishida, M.; Nakajima, T.; Honda, Y.; Kitao, O.; Nakai, H.; Klene, M.; Li, X.; Knox, J. E.; Hratchin, H. P.; Cross, J. B.; Adamo, C.; Jaramillo, J.; Gomperts, R.; Stratmann, R. E.; Yazyev, O.; Austin, A. J.; Cammi, R.; Pomelli, C.; Ochterski, J. W.; Ayala, P. Y.; Morokuma, K.; Voth, G. A.; Salvador, P.; Dannenberg, J. J.; Zakrzewski, V. G.; Dapprich, S.; Daniels, A. D.; Strain, M. C.; Farkas, O.; Malik, D. K.; Rabuck, A. D.; Raghavachari, K.; Foresman, J. B.; Ortiz, J. V.; Cui, Q.; Baboul, A. G.; Clifford, S.; Cioslowski, J.; Stefanov, B. B.; Liu, G.; Liashenko, A.; Piskorz, P.; Komaromi, I.; Martin, R. L.; Fox, D. J.; Keith, T.; Al-Lahan, A.; Peng, C. Y.; Nanayakkara, A.; Challacombe, M.; Gill, P. M. W.; Johnson, B.; Chen, W.; Wong, M. W.; Gonzalez, C.; Pople, J. A. Gaussian 03, revision C.02; Gaussian, Inc.: Wallingford, CT, 2004. (49) Becke, A. D. Phys. ReV. A 1988, 38, 3098. (50) Becke, A. D. J. Chem. Phys. 1996, 104, 1040. (51) Becke, A. D. J. Chem. Phys. 1993, 98, 5648. (52) Lee, C.; Yang, W.; Parr, R. G. Phys. ReV. B 1988, 37, 785. (53) Stephens, P. J.; Devlin, P. J.; Chabalowski, C. F.; Frisch, M. J. J. Phys. Chem. 1994, 98, 11623.

Sousa et al. (54) Wilson, P. J.; Bradley, T. J.; Tozer, D. J. J. Chem. Phys. 2001, 115, 9233. (55) Perdew, J. P. Phys. ReV. B 1986, 33, 8822. (56) Perdew, J. P. In Electronic Structure of Solids ’91; Ziesche, P., Eschig, H., Eds.; Akademie Verlag: Berlin, 1991; p 11. (57) Hamprecht, F. A.; Cohen, A. J.; Tozer, D. J.; Handy, N. C. J. Chem. Phys. 1998, 109, 6264. (58) Dunning, T. H.; Hay, P. J. In Modern Theoretical Chemistry; Schaefer, H. F., III, Ed.; Plenum: New York, 1976; pp 1-28. (59) Riley, K. E.; Op’t Holt, B. T.; Merz, K. M., Jr. J. Chem. Theory Comput. 2007, 3, 407. (60) Meijers, R.; Morris, R. J.; Adolph, H. W.; Merli, A.; Lamzin, V. S.; Cedergren-Zeppezauer, E. S. J. Biol. Chem. 2001, 276, 9316. (61) English, A. C.; Done, S. H.; Caves, L. S.; Groom, C. R.; Hubbard, R. E. Proteins 1999, 37, 628. (62) Kilshtain-Vardi, A.; Glick, M.; Greenblatt, H. M.; Goldblum, A.; Shoham, G. Acta Crystallogr., Sect. D: Biol. Crystallogr. 2003, 59, 323. (63) Hakansson, K.; Carlsson, M.; Svensson, L. A.; Liljas, A. J. Mol. Biol. 1992, 227, 1192. (64) Korndorfer, I. P.; Fessner, W. D.; Matthews, B. W. J. Mol. Biol. 2000, 300, 917. (65) Sankaranarayanan, R.; Dock-Bregeon, A. C.; Rees, B.; Bovee, M.; Caillet, J.; Romby, P.; Francklyn, C. S.; Moras, D. Nat. Struct. Biol. 2000, 7, 461. (66) Hall, T. M.; Porter, J. A.; Beachy, P. A.; Leahy, D. J. Nature 1995, 378, 212. (67) Cruickshank, D. W. J. Acta Crystallogr., Sect. D: Biol. Crystallogr. 1999, 55, 583. (68) Tobin, D. A.; Pickett, J. S.; Hartman, H. L.; Fierke, C. A.; PennerHahn, J. E. J. Am. Chem. Soc. 2003, 125, 9962. (69) Ryde, U.; Nilsson, K. J. Am. Chem. Soc. 2003, 125, 14232. (70) Schafer, A.; Horn, H.; Ahlrichs, R. J. Chem. Phys. 1992, 97, 2571. (71) Schafer, A.; Huber, C.; Ahlrichs, R. J. Chem. Phys. 1994, 100, 5829. (72) Alberico, W. M.; Garbarino, G. Phys. Rep. 2002, 369, 177. (73) Czuchaj, E.; Nicki, M. Chem. Phys. Lett. 2001, 335, 440. (74) Adamo, C.; Barone, V. Chem. Phys. Lett. 1997, 274, 242. (75) Koperski, J.; Czajkowski, M. Phys. ReV. A 2000, 62, Art 012505. (76) Wang, N. X.; Wilson, A. K. J. Chem. Phys. 2004, 121, 7632. (77) Zhao, Y.; Truhlar, D. G. J. Phys. Chem. A 2004, 108, 6908. (78) Riley, K. E.; Brothers, E. N.; Ayers, K. B.; Merz, K. M. J. Chem. Theory Comput. 2005, 1, 546.