CoâC Dissociation of Adenosylcobalamin (Coenzyme B12): Role of

Article pubs.acs.org/JPCA

Co−C Dissociation of Adenosylcobalamin (Coenzyme B12): Role of Dispersion, Induction Effects, Solvent Polarity, and Relativistic and Thermal Corrections Kasper P. Kepp* DTU Chemistry, Technical University of Denmark, Building 206, Kgs. Lyngby, DK-2800, Denmark S Supporting Information *

ABSTRACT: Quantum-chemical cluster modeling is challenged in the limit of large, soft systems by the effects of dispersion and solvent, and well as other physical interactions. Adenosylcobalamin (AdoCbl, coenzyme B12), as one of the most complex cofactors in life, constitutes such a challenge. The cleavage of its unique organometallic Co−C bond has inspired multiple studies of this cofactor. This paper reports the fully relaxed potential energy surface of Co−C cleavage of AdoCbl, including for the first time all side-chain interactions with the dissociating Ado group. Various methods and corrections for dispersion, relativistic effects, solvent polarity, basis set superposition error, and thermal and vibrational effects were investigated, totaling more than 550 single-point energies for the large model. The results show immense variability depending on method, including solvation, functional type, and dispersion, challenging the conceived accuracy of methods used for such systems. In particular, B3LYP-D3 seems to severely underestimate the Co−C bond strength, consistent with previous results, and BP86 remains accurate for cobalamins when dispersion interactions are accounted for.

■

± 5 kJ/mol in ethylene glycol30) and an activation enthalpy of ∼140 kJ/mol.31−33 However, the enzymes require reversible cleavage of this Co−C bond, a challenge solved by the molecular evolution of special properties of the enzymes and of AdoCbl itself.6,16 Several AdoCbl-dependent enzymes display catalytic rates (kcat) of ∼100 s−1, and Co−C bond cleavage is not rate-limiting since the enzymes shift the equilibrium constant between the Co(III) and Co(II) states by ∼1012, making it close to unity.31,34 This makes the reaction reversible and facilitates re-formation of the Co(III) resting state to complete the catalytic cycle. Consistent with this labilization, crystal structures that include substrates generally display broken Co−C bonds.35,36 This ingenious, controlled, and reversible production of an organic radical has been widely studied.3,5,6 It was previously thought that structural changes in the corrin could labilize the Co−C bond by steric strain, but this view has been largely discredited.6 The effect of the axial N-ligand is maximally 18 kJ/mol.37 Mutation of either the histidine that coordinates to cobalt in most enzymes or its hydrogen-bond partner Asp decreases turnover by only 1000-fold,38 a small fraction of the total catalytic effect and comparable to the effect of this ligand in the free cobalamins.37 Thus, the role of the Co−Nax bond is quite small. However, specific interactions between corrin and its side chains,31 tilting or pulling the Ado

INTRODUCTION Coenzyme B12 (5′-deoxyadenosylcobalamin, AdoCbl) is the largest and most complex cofactor in biology. It contains a central cobalt ion bonded to a true saturated carbon, rendering it the first discovered organometallic cofactor in nature.1,2 The unique ability of this cofactor to reversibly produce organic radicals for initiation of 1,2-shifts in substrates by homolytic Co−C bond cleavage has attracted considerable interest3−5 and, relevant to this work, has been a subject of intense theoretical investigation.6−22 AdoCbl includes a tetradentate corrin ligand binding to cobalt (Figure 1). Methyl, acetamide, and propionamide groups are linked to the corrin ring. One of the groups ends with a 5,6dimethylbenzimidazole (DMB) group (the nucleotide tail). AdoCbl is used as a cofactor by enzymes such as methylmalonyl coenzyme A mutase (MCAM),23 glutamate mutase (GLUMUT),24 class II ribonucleotide reductase, ethanolamine ammonia lyase,25 and diol- and glycerol dehydratase.26 It also serves a broad and vital role in the regulation of riboswitches.27 The enzymatic resting state contains Co(III). However, hemolytic cleavage of the Co−C bond of AdoCbl produces a five-coordinate Co(II) intermediate and an Ado radical that initiates radical reactions (1,2-shifts at carbon centers of the substrates).28 During the entire catalytic cycle, cobalt remains in the low-spin state.6 Most commonly (e.g., in human MCAM), a histidine coordinates to cobalt, producing the biochemically relevant imidazole-coordinated cofactor.23 AdoCbl is stable in solution due to a moderate Co−C bond dissociation enthalpy (BDH) (∼126 ± 8 kJ/mol in water,29 132 © 2014 American Chemical Society

Received: April 13, 2014 Revised: June 23, 2014 Published: August 12, 2014 7104

dx.doi.org/10.1021/jp503607k | J. Phys. Chem. A 2014, 118, 7104−7117

The Journal of Physical Chemistry A

Article

Figure 1. Studied model of histidine-bound AdoCbl, optimized in a Cosmo model. More details, including full coordinates, are given in the Supporting Information.

group,39 and electrostatic favoring of the dissociated state17,22 may affect Co−C bond cleavage. Since the side chains in AdoCbl have never been modeled in total quantummechanically, their effect on Co−C bond cleavage remains indirectly determined from, e.g., QM/MM partition calculations.17 To understand the nature of the Co−C bond and its cleavage, density functional theory (DFT)40,41 has been widely applied since the beginning of the millennium7,8,12,44 to study corrin models.42−54 The first of these studies used the B3LYP functional,55,56 which is the most widely used functional today57 and has worked successfully for other systems.58,59 However, B3LYP was found to underestimate Co−C bond dissociation energies (BDEs) by 35−80 kJ/mol,42 a potential major problem if common to other systems, since errors in metal− ligand bonds will directly influence reactions where such bonds form or break in catalysis and bioinorganic chemistry. Investigation of metal−carbon BDEs in related systems15 indicated that the Hartree−Fock (HF) exchange in B3LYP and, to some extent, the LYP correlation functional56 were responsible for the small BDEs, and the non-hybrid BP8660,61 has since been the preferred method for theoretical B12 chemistry.15 Subsequently, benchmarking of B3LYP was done against available experimental BDEs for M−L diatomics where dispersion and solvent effects are conveniently negligible, for both the first row62,63 and the second and third rows of the dblock.64 Across the d-block, B3LYP under-binds and BP86 over-binds by similar amounts for the diatomics, and HF exchange could largely explain the reduction in BDEs of metal− ligand bonds, suggesting the use of 10% hybrids such as TPSSh,65,66 which was later confirmed to be most accurate across the d-block for describing M−L diatomic bond strengths.63,64 To understand the Co−C cleavage process, in particular the role of the side chains, and to evaluate in this light the performance of B3LYP and BP86, this paper reports computations of Co−C cleavage including all the side chains contributing to the dispersion effect on Co−C cleavage, leaving out only the large, downward-pointing nucleotide tail. Physical and chemical effects on Co−C cleavage were investigated, giving the most complete view of this process in the real coenzyme. Dispersion−solvent interactions substantially shape the potential energy surface (PES) of Co−C bond cleavage. The study finds that BP86-D3 is more accurate for modeling

Co−C cleavage than B3LYP-D3. In contrast to the studies of diatomics,62,63 AdoCbl requires less than 10% HF exchange, probably due to its closed-shell bond breaking, which is not typical across the d-block M−L diatomics. The results confirm previous findings from M−L diatomics62 and porphyrins,67−69 showing that B3LYP under-binds substantially even with dispersion corrections. However, the study also details how large, realistic cluster models suffer from large, compensatory dispersion−solvent interactions. Since the calculations include all relevant interactions with the dissociating group, this nonenzymatic reference system thus provides a good starting point for future studies of enzymatic Co−C cleavage.

■

METHODS Geometry Optimization. All calculations were performed with the Turbomole 6.3 software.70 The crystal structure of AdoCbl has the nucleotide tail with its terminal DMB coordinating to cobalt.71,72 In the calculations, the nucleotide tail was changed to only extend to the amide functional group, as this large group points downward, away from the Ado group that is cleaved and has little role in the Co−C cleavage both with and without enzymes.37 The model had a total charge of +1 and included 175 atoms. Optimization of all geometries were performed at the BP86/ def2-SVP level73 using the Cosmo solvation model74 since condense phase screening has been shown to improve the accuracy of geometries of charged systems vs experimental data.75 The default optimized radii of 2.0 Å for C, 1.83 Å for N, 1.72 Å for O, 1.3 Å for H, and 2.0 Å for Co. (Sensitivity to this choice was investigated separately, as described below.) Structural sensitivity to electrostatic screening was investigated by performing geometry optimizations of all species (complex and dissociated states) with dielectric constants of 1, 10, and 80, since the geometry of a large molecule with many soft interactions might be sensitive to changes in solvent polarity. Electronic energies were converged down to 10−7 hartrees, and to achieve stable equilibrium geometries, the gradient was converged down to 10−3 a.u. Calculations of Co−C BDEs. The Co−C bond strength was computed using eq 1: BDE = E[Co(II)] + E[Ado•] − E[Co(III)Ado] + ΔZPE + ΔHcorr 7105

(1)



Article

multiple interactions between Ado and corrin, and constitutes a prominent example of the full range of challenges associated with modeling large inorganic systems. To the computed electronic energies were added corrections for scalar-relativistic effects calculated with the Moloch module of Turbomole. Scalar-relativistic effects may amount to ∼10 kJ/ mol or more in first-row transition metals and can, in some cases, be systematic.57 Energies were computed both with and without D3 dispersion corrections, to evaluate their functionalspecific effects. As with the geometry optimization, the singlepoint energies used for computing the BDEs were evaluated using dielectric constants of 1, 10, and 80, and by varying the radius of cobalt, because this atom is exposed during Co−C cleavage; thus, it was thought that solvent probing of this exposure might contribute to solvent effects on the computed BDE. Choice of Functionals. Computational chemistry investigations should routinely use several density functionals to probe sensitivity of results to functional type, since variable HF exchange and other ingredients of the functionals substantially affect the computed relative energies.57 In this work, the BDEs were computed with functionals B3LYP,55,56 BP86,60,61 BLYP,56,60 PBE,83 PBE0,83 and TPSSh,65,66 using the def2TZVPP basis set,84 and were converged to 10−7 hartrees. The chosen functionals represent a range of HF exchange from 0% (BP86, BLYP, PBE) to 10% (TPSSh), 20% (B3LYP), and 25% (PBE0). As a good practice, B3LYP should always be included simply due to its order-of-magnitude dominance in the field; i.e., it is the current standard despite known shortcomings.57 BP86 should be included because it has been the standard of cobalamin chemistry since 2003,15 providing the most accurate BDEs and geometries in many studies.6,10,17,45,53,54,85,86 Furthermore, comparing BP86 vs BLYP allowed the probing of the role of, i.e., the LYP correlation functional previously found to lower BDEs of organometallic bonds.15 As an additional feature, PBE, PBE0, and TPSSh were parametrized using fewer parameters and more fundamental constraints.83,65,66 TPSSh was included for its reported accuracy in other first-row transition-metal systems,57,87,88 which was further documented for spin crossover of cobalt and iron systems67,89 and metal−ligand bond dissociations63,64 relevant to the present study. Calculation of Co−C Bond Dissociation Potential Energy Surfaces. The full PESs were also computed for the Co−C process, to understand the cleavage process in detail. The PESs were computed by geometry optimizing the electronic configuration with the Co−C distance fixed at 20 points in the range from 1.9 to 5 Å. The Co−C distance is 3.66 Å in the human enzyme MCAM which is also histidine-on, nucleotide tail off, based on structure 2XIQ.90 As for the geometry optimizations used for computing the formal BDEs, the method was BP86/def2-SVP. The electronic energies in hartrees and relative energies in kJ/mol are given in Supporting Information, Tables S1 and S2, to enable reproducibility, and the BP86-optimized PES is given in Supporting Information, Figure S1. Subsequently, at each geometry-optimized position of the PES, single-point def2-TZVPP energies were computed with B3LYP, TPSSh, BP86 with dielectric constants of 1, 10, or 80, with or without dispersion corrections (a total of 12 types of calculations and 240 converged def2-TZVPP single-point energies for the PESs). Since the molecular orbitals come closer in energy as the bond length is extended, all

Here, BDE is the homolytic bond dissociation energy, E[Co(II)] is the computed quantum-mechanical energy of the fully geometry-optimized Cob(II)alamin model, E[Ado•] is the energy of the isolated adenosyl radical, E[Co(III)Ado] is the energy of the Co(III) AdoCbl model, ΔZPE is the differential vibrational zero-point energy, which favors dissociation,57 and ΔHcorr is a correction term of thermal and enthalpy effects, calculated from the vibrational state function for all three species. No scale factors were applied as BP86 has optimal vibration scale factors very close to unity (0.996) for the studied type of fundamental vibrations.76 In total, 36 models were used to compute the energy at each of the three optimized geometries for each species, giving 324 converged DFT/def2-TZVPP single-point calculations for this part of the work (4538 basis functions for AdoCbl model, 3270 for Co(II)Cbl model). Counterpoise corrections to the basis set superposition errors (BSSEs) were computed using the jobbsse script of Turbomole, using the optimized AdoCbl structure divided into Ado (30 atoms) and resulting Cob(II)alamin (145 atoms) fragments for the calculation (Cob(II)alamin: atoms 1− 71, 90−162, and 175; see Supporting Information for xyz coordinates of the full model). Many transition-metal systems are open-shell before metal− ligand bond dissociation, reducing the effect of HF exchange, explaining why B3LYP has not been highly problematic throughout inorganic chemistry.57 However, the Co(III) state of cobalamins is closed-shell, and therefore the HF bias toward the dissociated state is larger than in other systems.15 Furthermore, the size and complexity of cobalamins produce large dispersion and solvent effects that complicate further their modeling. Accordingly, the finding that BP86 outperforms B3LYP15 for cobalamins has recently been challenged.20,21 This challenge mainly finds its support in the advent of dispersion corrections to DFT. In particular, DFT-D3 by Grimme and coworkers77 has been applied with good experience.57,78 Dispersion increases BDEs of bulky molecules,57 as also seen by recent computations of the full methylcobalamin (MeCbl).79 Therefore, it was suggested that the accuracy of BP86 is due to its over-binding tendency mimicking dispersion (cancellation of errors).20,21 Using the strongly binding early dispersion correction of Grimme80,81), B3LYP becomes more accurate.21 However, D3, the recent, recommended dispersion correction, is markedly smaller for bulky transition-metal systems than D2.18 Thus, in other reports, BP86 remains accurate when dispersion is included for similar models, whereas B3LYP underestimates the Co−C BDE substantially.82 B3LYP-D2 produced 36 kJ/mol higher BDE in MeCbl than B3LYP-D3.79 Thus, whether dispersion “saves” B3LYP when describing metal−ligand bonds such as Co−C remains unclear. Furthermore, the diatomics with negligible dispersion effects suggest that 10% is most accurate;63,64 i.e., an accurate dispersion model such as D3 together with 10% HF exchange TPSSh could give the best results. Altogether, there are still varying opinions, and the cobalamins are currently at the center of some of these debates.18,20,21 AdoCbl has stronger dispersion interactions than MeCbl due to its large Ado group: Previous work found B3LYP*-D2 to predict Co−C BDEs accurately.20 However, as other previous work, the used AdoCbl models neglected side chains. Thus, agreement with experimental data where these side chains are present could be due to a cancellation of errors from the too strongly binding early dispersion corrections80,81 and absence of corrin side chains. The real AdoCbl cofactor is crowded with 7106



Article

being the softest in the first coordination sphere and many times softer than the Co−C bond, and thus more sensitive to environmental changes.92 These minor geometric changes in the first coordination sphere are not enough to explain the substantial vacuum/solvent effect on Co−C cleavage observed experimentally for cobinamides.54 In Table 1, the computed structures of the AdoCbl and cob(II)alamin models have been compared with experimental structures of AdoCbl71,72 and cob(II)alamin.93 One can deduce a mean absolute error for all metal−ligand bonds of 0.03 Å for both AdoCbl and Cob(II)alamin. The computed Co−C bond length, which is critical to the present work, is within 0.02 Å of available experimental structural data. An important change in geometry between the Co(III) and Co(II) states is the experimental 0.01 Å shortening of the average equatorial Co−Neq bonds in Cob(II)alamin relative to AdoCbl. Taking the average of the four bonds in the three computed structures, this bond contraction is modeled to be almost exactly 0.01 Å, showing the strength of computational chemistry for relative trends rather than absolute magnitudes. This contraction may affect the Co−C BDEs and is thus important to reproduce. The softer Co−Nax bond contracts by 0.1 Å upon Co−C cleavage according to the crystal structures; the computed structures suggest a contraction of 0.05 Å, the difference being large because the bond is soft and hence more variable: This is also seen from the Co−Nax bond length from neutron diffraction data of AdoCbl at 15 K,71 which is close to the computed “0 K” bond length, with only a 0.02 Å difference. Several amide side chains interact with the dissociating Ado group, as seen in the optimized structure in Figure 1. Two of them form direct hydrogen bonds with Ado. These hydrogen bonds, present in the real cofactor, will increase the BDE relative to simpler corrin models. In a real, condense phase medium, solvent hydrogen bonds will compensate the loss of the hydrogen bonds. Also, one of the methyl groups that points downward (toward the α-side, or N-base side, of the cofactor) is very close to cobalt (hence it is referred to here as the proximal methyl group) and modulates the conformation of the corrin ring and the axial N-base, as seen in the crystal structure90 and also in the BP86-Cosmo optimized structure: The proximal methyl group causes the axial imidazole to tilt slightly, an effect that, despite the crystal structure evidence, has not been discussed widely and also could impair the accuracy of truncated corrin models. As shown in Figure 1 and in detail in the Supporting Information, the calculations reproduce such a tilt, although it again does not impact the Co−C bond. Co−C Bond Dissociation Energies: Effects of Functional, Dispersion, and Solvent. After validating the properties and accuracy of the computed structures, the BDEs as computed directly on these structures in the three types of environment (ε = 1, 10, and 80) were compared. The functionals BP86, BLYP, B3LYP, PBE, PBE0, and TPSSh were used both with and without dispersion and with all three values of ε, to compute the energies of AdoCbl, cob(II)alamin, and Ado structures and calculate the BDEs. These results are reported in Figure 2 without any corrections for ZPE, thermal effects, and scalar relativistic effects, in order to first address the effects of functional, dispersion, and solvent directly (numerical data can be found in Supporting Information, Table S18). Before discussing the results, the corresponding BDEs calculated from truncated corrin models allow a quantification of the importance of the side chains in AdoCbl. The Co−C BDE from truncated corrin with Ado and imidazole as axial

computations were first done for the triplet states, and then the open-shell singlet states that correspond to the brokensymmetry solutions to the true singlet state were converged starting from the triplet-state orbitals to avoid local minima, i.e., the closed-shell singlet. Scalar-relativistic effects, which can be of the order of ∼10 kJ/mol for M−L bond energies of first-row transition metals,57 have been computed for all configurations along the binding curves, using the Cowan−Griffin operator91 as implemented in Turbomole. Additional numerical data can be found in the Supporting Information relating to the effect of Cosmo probe radius (Table S3), vibrational and thermal corrections (Tables S5 and S6), energies in hartrees for the BDE calculations (Tables S7−S9), and relativistic corrections as output hartrees (Table S10). A structural alignment of the optimized Co(II) and Co(III) states is given in Figure S2, and the xyz coordinates of optimized structures are given in Tables S11−S16. Details on solvent effects are given in Tables S17 and S18.

■

RESULTS Optimized Geometries of AdoCbl and Cob(II)alamin. The molecular geometry of the full model of AdoCbl was optimized both in vacuum and using Cosmo models of condense phase since this was anticipated to affect the structure of a large bulky molecular system as AdoCbl. The optimized geometry in a water-like solvent model is shown in Figure 1; the full xyz coordinates are given in the Supporting Information, Tables S11−S16. The solvent affects the sidechain geometries and the details of the Ado orientation, but the first coordination sphere is essentially unaffected by the solvent, as shown in Table 1. The only significant change occurs in the Co−Nax bond, which is elongated by 0.01−0.02 Å when optimized in water. This change is consistent with this bond Table 1. Computed vs Experimental Metal−Ligand Equilibrium Bond Lengths (Å)a this work previous work

X-ray

neutron diffraction

2.024b 1.859b 1.859b 1.932b 1.932b 2.155b

2.03c 1.88c 1.87c 1.91c 1.92c 2.24c

2.02d 1.89d 1.89d 1.91d 1.91d 2.21d

Co(II)Cbl 1.878 1.901e 1.883 1.936 1.951 2.161 2.154e

1.87f 1.87f 1.89f 1.91f 2.13f

ε=1

ε = 10

ε = 80

Co−C Co−N1 Co−N2 Co−N3 Co−N4 Co−Nax

2.035 1.889 1.908 1.947 1.952 2.190

2.032 1.888 1.907 1.946 1.953 2.200

AdoCbl 2.032 1.888 1.906 1.945 1.954 2.206

Co−N1 Co−N2 Co−N3 Co−N4 Co−Nax

1.879 1.884 1.936 1.952 2.154

1.878 1.884 1.937 1.951 2.160

a

For system comparison, all computational data reported were done with BP86. ε refers to the dielectric constant of the medium used for optimizing the geometry of the AdoCbl model. bDMB-on model of AdoCbl from ref 86. Co−N1 and Co−N2, and Co−N3 and Co−N4 are averaged. cFrom ref 72. dFrom ref 71. eTruncated corrin model of cob(II)alamin from ref 15. fAverage Co−Neq distances and Co−Nax From ref 93. N1 is the Neq closest to proximal methyl (N21 in ref 93.); N2 is the second-closest to proximal methyl (N24 in ref 93.), followed in the same direction by N3 and N4. 7107



Article

Figure 2. Computed Co−C bond dissociation energies of AdoCbl on structures optimized with three different dielectric constants.

dispersion is even more important, and although D3 corrections are smaller than the earlier dispersion corrections,18 they still exceed 100 kJ/mol when all the relevant interactions in AdoCbl are included. The conclusion is thus that DFT without dispersion is meaningless for studying bond breaking and formation reactions of large crowded systems such as cobalamins, but another major issue is that these effects are not balanced by explicit solvent dispersion terms.18 This problem is one of the most severe ones facing current quantum-chemical cluster modeling. (ii) Except for a very few cases, the structures optimized in water consistently give larger BDEs than structures optimized in vacuum, with the differences being largest for the singlepoint vacuum energies; i.e., the geometric differences are reduced once a Cosmo model is applied also to the single-point energies. However, this effect is relatively minor and will not be discussed further in the following. (iii) The solvent effect on the computed BDE is very large, typically 70−80 kJ/mol, and much larger than solvent effects found from previous truncated corrin models. The most important lowering in BDE comes moving from ε = 1 to 10, consistent with the most electrostatic screening effect occurring at small values of the dielectric constant, since the interactions scale inversely with ε. These large solvent effects are consistent with the solvation penalty of hydrating the Ado radical, as also seen experimentally, both when comparing AdoCbl in water and ethylene glycol (6 kJ/mol),29,30 and in particular, when comparing cobinamides in gas phase and water (30 kJ/mol).54 The nature of the solvent is also important, with a computed reduction in BDE from ε = 10 to 80 of 10−12 kJ/mol. The chemical variation between cobalamins and cobinamides is unlikely to explain the difference. Rather, both the approximations of the Cosmo model and the fact that the gas phase experiment may have some shielding effects from the applied apparatus and impurities could contribute. Also, the experimental water solutions contained phosphate buffer.37 H2PO4−

ligands was computed with BP86 (no dispersion) to be 142.5 kJ/mol, i.e., smaller than with only the ribose ring (155 kJ/mol) or methyl (160 kJ/mol).17 This showed that there is a considerable steric repulsion between the Ado group and the corrin, which outweighs the increased bonding caused by hydrogen bonds. However, in addition to not including quantum-mechanically the side chains, these numbers were also without dispersion, which should substantially influence the BDE of crowded systems like AdoCbl.18,20,21 From the computed BDEs in Figure 2, several observations can be made: (i) For the large AdoCbl model, any density functional that does not include dispersion fails substantially due to underbinding. This finding is fully in line with recent benchmarks on methylcorrin using a range of functionals, where only dispersion-corrected methods succeeded in producing BDEs resembling renormalized coupled-cluster and experiment.82 However, the under-binding of many functionals is even more dramatic than previously observed because side chains interacting with Ado have not been included quantummechanically before. Even without dispersion, B3LYP predicts binding for a truncated corrin model of AdoCbl,20 but in the Cosmo models of the large model used here, many of the functionals do not predict binding; i.e., the BDE is negative. As the side chains have repulsive parts included in the functional, but not the attractive van der Waals interactions, a very large and crowded system such as AdoCbl is meaningless to describe without dispersion. As an example, BP86 without dispersion in vacuum predicts a BDE of 71−79 kJ/mol, but the value in water is only 3−4 kJ/mol, with the experimental BDH in water being 126 kJ/mol.29 The enthalpy corrections to the energies are of the order of ∼12−13 kJ/mol (see below), so they cannot explain this enormous difference. For truncated models studied before, the D2 dispersion correction gave ∼50 kJ/mol larger BDEs, which ended up giving reasonable results for B3LYP.20 For our large model, 7108



Article

Because of the size of the system studied (175 atoms) and the close interactions between several side chains and the Ado group, the total complex benefits substantially from the presence of basis functions from other atom centers. This gives rise to non-negligible corrections, despite the use of the def2-TZVPP basis set. Whereas for truncated simple methylcorrins, this correction amounts to 3−6 kJ/mol with polarized triple-ζ basis sets,15,82 it ranges from 9 to 12 kJ/mol in the AdoCbl model with side chains present, depending on the nature of the functional used, showing the important role of both number of atoms and basis functions. In fact, the BSSE seems to grow roughly proportional to N/b, where N is the number of interacting atoms and b is the number of basis functions. The three hybrid functionals provide the smallest corrections, consistent with their tendency of producing more diffuse electron densities;57 thus, hybrids are less prone to BSSE because their electron densities are already more polarized by the HF exchange of the functional itself, thus gaining less from the spatial freedom of additional ghost orbitals. The BSSE is very system-size dependent and is not generally included during parametrization of standard functionals to experimental thermodynamic data, which is a main reason why counterpoise corrections often impair results. (This argument is more valid for DFT than for non-empirical correlation methods, viz., the “over-correction debate”.96) Thus, standard research papers rarely report corrected energies for transition states and chemical energies in general. When included, the bound state is disfavored by additionally 10 kJ/mol to further increases the preference toward non-hybrid dispersion corrected methods. Relativistic, Thermal, and ZPE Corrections to the Co− C BDEs. To judge the performance of functionals effectively, one should correct for the chemical and physical context of the experimental measurements, which is not always done. For example, the accuracy of functionals for estimating the ground state spin state of a complex should include the vibrational entropy which favors high-spin by typically ∼25 kJ/mol at room temperature; i.e., accurate functionals should actually give energies that favor low-spin by a similar amount.57 Thus, lack of attention to important systematic thermodynamic, zero-point energy, or relativistic effects leads to erroneous benchmarking results. Table 3 shows the computed BDEs and BDHs of the six applied DFT-D3 methods corrected for ZPE, thermal and scalar-relativistic corrections. Consistent with Figure 2, the nonrelativistic BDEs in water vary from 63 kJ/mol with PBE0-D3 (25% HF exchange) to 140 kJ/mol with BP86-D3 (0% HF exchange). The scalar-relativistic corrections to the BDEs are relatively minor in these cases but grow rapidly down the periodic table: they amount to 3−6 kJ/mol and strengthen the Co−C bond. For truncated methyl corrins, the relativistic corrections increase BDEs by 6−7 kJ/mol;15 i.e., large and small models give similar relativistic corrections, since the corrections mainly involve the core electrons on cobalt during Co−C cleavage, which is less affected by substituents and the nature of the dissociating alkyl group. As explained previously,57 scalarrelativistic effects on M−L bonds follow the simple principles from atomic orbitals;98 i.e., they strengthen σ-bonds due to relativistic contraction, whereas higher angular momentum molecular orbitals in bonds with bond orders >1 are destabilized and expanded.57

is a strong salting-out Hofmeister ion that reduces the apparent solubility of solutes and polarity of water.94 Thus, true water may have an even smaller BDE than that reported, but this remains to be investigated. (iv) The last observation from Figure 2 relates to how the functionals perform. To address this question accurately, one needs to include additional corrections, i.e.. ZPE, thermal corrections to the enthalpy, and scalar-relativistic corrections, that all affect computed BDEs.15 However, before correcting for these effects, some trends are already deducible from Figure 2. The most important is that all functionals, once corrected for dispersion and solvent effects, produce moderate binding energies consistent with an organometallic bond (50−200 kJ/ mol). The trend is PBE0 ∼ B3LYP < TPSSh < BLYP < PBE < BP86. The PBE0 and B3LYP functionals produce the smallest BDEs, consistent with their larger components of HF exchange (25% and 20%) which polarizes metal−ligand bonds and reduces bond strengths.62 It can be seen that the LYP functional, which is known to reduce bond organometallic strengths,15 works to the effect of ∼5% HF exchange, bringing B3LYP and PBE0 very close in performance. Then follows TPSSh (10% HF exchange), BLYP (0% HF exchange but with the BDE-reducing LYP functional, confirming that this effect resembles 5% HF exchange as BLYP lies equidistant between BP86 and TPSSh), PBE, and BP86. Except for the subtle differences between last two nonhybrid functionals, these trends are well understood and important to consider in computational inorganic chemistry.57 The observation that B3LYP performs poorly even with dispersion corrections is consistent with recent benchmarks of methylcorrin BDE.82 Given the order-of-magnitude dominance of B3LYP as the current standard in the field,57 these works together show that for organometallic bonds resembling Co−C, B3LYP, even with dispersion corrections, should be avoided. Correction for Basis Set Superposition Error. An important contribution to errors in computational chemistry is the use of limited basis sets. In calculation of interaction energies such as BDEs, the electrons of the complex benefit from basis functions centered on atoms that belong to both fragments, whereas the fragments do not contain these additional functions.95 This causes an additional variational freedom and resulting stabilization of the electrons in the complex, leading to a bias toward the complex known as the BSSE. To estimate the contribution of this error to the computed BDEs, counterpoise corrections96 according to the method of Boys and Bernardi97 were computed using the jobbsse script of Turbomole. Since the BSSE might depend on the functional used, corrections were computed for all six functionals used in this work. The corresponding results are shown in Table 2. Table 2. Computed Counterpoise Corrections of the Basis Set Superposition Error

functional

HF (%)

corrected electronic energy (a.u.)

BSSE correction (a.u.)

effect on computed BDH (kJ/mol)

PBE0 B3LYP TPSSh BLYP PBE BP86

25 20 10 0 0 0

−5455.4611 −5457.8575 −5460.5657 −5458.9049 −5455.2650 −5460.8796

0.0033 0.0036 0.0033 0.0044 0.0043 0.0037

−8.7 −9.5 −8.6 −11.6 −11.2 −9.7 7109



Article

Table 3. Computed Co−C Bond Dissociation Energies (BDEs) and Enthalpies (BDHs) without or with Scalar-Relativistic Corrections (in kJ/mol) functional

HF (%)

BDEnon‑rel (ε = 80)

BDErel (ε = 80)

rel corr

BDHrel (ε = 80)

BDHrel (ε = 41)

comp − expta

PBE0-D3 B3LYP-D3 TPSSh-D3 BLYP-D3 PBE-D3 BP86-D3

25 20 10 0 0 0

62.6 65.9 89.7 110.2 119.8 139.6

68.4 71.0 94.3 113.7 123.3 143.0

5.8 5.1 4.6 3.5 3.5 3.5

55.8 58.4 81.7 101.1 110.7 130.4

57.6 60.2 83.5 102.9 112.5 132.2

−70.2 −67.6 −44.3 −24.9 −15.3 4.4

a Computed minus experimental BDH, including all corrections, in water (ref 29). These numbers are 4 kJ/mol more negative for ethylene glycol (ε = 41). The experimental error bar is 8 kJ/mol (ref 29).

Figure 3. (A) Potential energy surfaces of Co−C cleavage in AdoCbl, computed with several methods. (B,C) Structural changes occurring at Co−C = 2.8 Å (B) and 3.5 Å (C).

In contrast to the relativistic corrections, which increase bond strengths, the ZPE and enthalpy corrections lower M−L bond strengths: the combined effect means that BDEs should be corrected by approximately 12−13 kJ/mol to yield the BDHs measured experimentally. When scalar-relativistic and thermalvibrational effects are both added, they cancel partly to yield

relatively minor total corrections to the naively computed BDEs. With these corrections, BDHrel, the scalar-relativistic corrected BDH including thermal and vibrational ZPE corrections in polar solvents with ε = 80 or 41, represents the relevant estimate for properly comparing the functionals, since any other choice would reflect cancellation of errors.57 7110



Article

Figure 4. Scalar-relativistic corrections to the Co−C cleavage process (REL = relativistic, N-REL = non-relativistic).

sponding to Co−C at 2.8 and 3.5 Å shown in Figure 3B,C (see numerical data in Supporting Information, Tables S1 and S2). From the PESs in Figure 3A, one can deduce several conclusions on the importance of dispersion and solvent, as well as functional consistent with those obtained for the computed BDEs: BP86 gives stronger binding than TPSSh, which gives stronger binding than B3LYP, and B3LYP without any dispersion does not predict any binding at all. With all side chains included that interact with Ado, the dispersion can be seen to be crucial in shaping the PESs. All PESs show a local minima ∼2.8 Å that relates to a particularly favorable orientation of the Ado group once it has dissociated, by orienting toward the amide groups of the cobalamin but still with incomplete radical formation, as seen from the dissociating carbon being between sp2- and sp3hybridized in its structure (Figure 3B; the radical formation is explored further below). From this minimum, the dispersioncorrected methods increase the energy as Co−C breaks, because the favorable short-range dispersion interactions are reduced as Ado moves beyond this minimum. This minimum is not likely to be experimentally observable because it is shallow relative to the equilibrium Co−C bond distance; i.e., it will have a short lifetime. Also, this is an enthalpy minimum that may be further blurred by thermal motion in a real condense phase. Thus, the experimental BDH and activation enthalpies refer to estimates of fully dissociated species, as appropriate. The difference between the energy of the PES at a minimum representing an observable species and the BDE calculated from isolated reactants and products has been referred to as a cage effect in enzyme models, estimated to be 20−25 kJ/mol for truncated corrins without dispersion.17,19 A cage is evident in the crystal structure of MCAM, corresponding to a Co−C distance of 3.66 Å (PDB: 2XIQ).90 The Ado radical is fully formed at 3.5 Å, as judged from the fully sp2-hybridized Ado

The last column of Table 3 reports this relevant difference between computed and experimental BDHs in water. In ethylene glycol, the corresponding numbers have 4 kJ/mol larger errors (except for BP86-D3, which becomes essentially zero) since the computed effect of moving from water to ethylene glycol is +2 kJ/mol, but the experimental effect is +6 kJ/mol,29,30 probably due to deficiencies in the implicit solvation models (see Supporting Information, Table S17). Thus, using this comparison of the six methods, BP86-D3 and to some extent PBE-D3 perform very well also for the large AdoCbl model, despite dispersion now being included for all interactions. Thus, the previously reported accuracy of BP86 for truncated corrins,15 suggested to be due to cancellation of overbinding and lack of dispersion,20 was most likely rather due to cancellation between side-chain neglect and the total molecular dispersion energy correction. As a final test of the model, since the solvent effects were found to be very large, the Cosmo radius of cobalt was varied, because it is exposed during Co−C bond dissociation and thus would contribute selectively to the solvation energy of the dissociated Co(II) state. However, a total variation of only 3 kJ/mol in BDE was caused by the choice of Co solvation radius in the realistic range from 1.5 to 2.5 Å (see Supporting Information, Table S3), something not previously investigated but potentially important when atoms are exposed in solvent during a computed reaction. Thus, in this case, this radius effect is minor and will not be discussed further in the paper. Co−C Potential Energy Surfaces of AdoCbl. The computed, full PESs of Co−C cleavage for the large model are given in Figure 3A. The computed PESs with BP86, TPSSh, and B3LYP both with and without the most recent D3 correction enable an evaluation of the magnitude of dispersion and side-chain interactions during the Co−C cleavage process of the non-enzymatic coenzyme. Relevant structures corre7111



Article

Figure 5. Spin density plot showing the radical formation upon Co−C cleavage of AdoCbl, computed with BP86-D3, TPSSh-D3, and B3LYP-D3 using the def2-TZVPP basis set.

relativistic effects tend to increase the energy of the fragments, and thus, the BDE, consistent with the findings obtained for the isolated species in Table 3. Radical Formation Is Complete at 3.5 Å. To add further details to the computations, the electronic structure upon Co− C cleavage was analyzed in more detail. The biologically critical function of AdoCbl is to produce reversibly the radical needed for initiating radical formation and subsequent 1,2-shifts within the enzymatic substrate.3,6 Enzymes such as GLUMUT and MCAM have characteristic EPR spectra that are interpreted as having the substrate-based radicals separated by 5−7 Å in the proteins, consistent with the Ado carbon located in between perhaps 4 Å from Co. The lifetime of the Ado radical is too short to be measured.4 A relevant question is at which Co−C distance the Ado radical is fully formed. This has been investigated in detail before by computing the spin polarization corresponding to the separation of the α and β electrons along the Co−C PES,17 giving a plot of the type shown in Figure 5. For the truncated corrin, the spin polarization started at 2.4 Å and was complete at 3.5 Å (BP86 functional).17 In Figure 5, the corresponding spin polarization has been compared now for the large AdoCbl model using BP86-D3, TPSSh-D3, and B3LYP-D3 with the def2-TZVPP basis set. It can be seen that the two hybrid functionals follow similar tendencies with radical formation beginning at 2.6 Å, 95% complete at 3.0 Å, and 100% complete at 3.5 Å. However, with the BP86-D3 method, the radical first emerges at 2.7 Å, builds up more slowly and is only 80% complete at 3.0 Å, but then is virtually complete at 3.5 Å. The fact that the hybrid functionals generate spin polarized electron structures at shorter Co−C distance is consistent with the known tendency of HF exchange in hybrids to polarize electron and spin densities:57,62 This is the same physics that increases bond lengths and reduces bond strengths. It is notable that the point of full radical formation, 3.5 Å, is insensitive to method and quite close to the observed

carbon (Figure 3C), indicating that this distance in the enzymes is enough to enable radical formation (confirmed below using spin densities). Thus, this calculated 3.5 Å value serves as a minimum threshold for the Co−C distance required by AdoCbl enzymes to enable radical formation and initiation of subsequent reactions. Consistent with this result, in all reported crystal structures, the Co−C distance is close to or larger than this distance. At this Co−C distance, BP86-D3 gives an energy of ∼70 kJ/ mol (∼63 kJ/mol with corrections) vs a corrected BDH of ∼130 kJ/mol (Table 3), a 67 kJ/mol difference. The deduced experimental shift in the equilibrium constant of Co−C cleavage toward unity amounts to a factor of 1012 or a ΔG ≈ 70 kJ/mol, with entropy being the minor part in the enzyme31,34 and TΔS ≈ 25 kJ/mol in MCAM.99 Thus, ΔH for Co−C cleavage is ∼25 kJ/mol as well, due to Keq ≈ 1. Using the computed 63 kJ/mol at 3.7 Å, the enzymes thus seem to favor the “cage” Ado radical state by ∼38 kJ/mol. These numbers (of course with sizable error bars) should be of relevance to future discussions on the enzymatic Co−C bond labilization. This energy estimate is not very sensitive to functional type (but is, of course, to dispersion): TPSSh-D3 and B3LYP-D3 with ε = 10 thus give 40 and 10 kJ/mol at 3.7 Å from Figure 3, vs 90 and 66 kJ/mol in water uncorrected, a difference of at least 50 kJ/mol. If dispersion is removed, BP86 (ε = 10) gives 40 kJ/mol. The large effects when including all dispersion interactions in AdoCbl is in contrast to previous work that did not treat side chains quantum-mechanically.14,17,22 Furthermore, as already partly evident from the relativistic corrections to the BDE calculations in Table 3, when computing scalarrelativistic corrections for the entire PES, the effects are minor, within 4−5 kJ/mol, and quite systematic, as seen in Figure 4. From the relativistic and non-relativistic PES, it can be seen that relativistic corrections only play out upon bond dissociation and mainly when the radical has been formed, and the 7112



Article

hybridized geometries and from spin densities. Thus, the experimentally observed enzyme geometries are quantummechanically sufficient to produce the Ado radical in the enzymes, which is then required in the catalytic reactions. These facts should be of relevance for understanding also how the enzymes labilize this enigmatic organometallic cbond. The choice of accurate models is crucial to catalytic processes where metal−ligand bonds form or break. An important challenge in the theoretical understanding of Co−C labilization is to separate the effect of the cofactor itself from that in the full enzyme. AdoCbl is large (175 atoms in the model without the nucleotide tail as used here) and has not been described fully quantum-mechanically before, despite the advent of dispersion corrections that render such calculations on large systems meaningful. In contrast, DFT-D calculations of the analog MeCbl, which has substantially less complicated interactions with the dissociating methyl radical, were recently reported.20,21 The present work indicates that the functional BP86, suggested previously to outperform other methods for this problem,15 remains the most accurate functional when the full interactions are included quantum-mechanically. The insight that B3LYP underbinds both from metal−ligand bond strengths in general62,63 and from recent large-scale benchmarks of the methylcorrin BDE,82 is confirmed in the present work where the most recent D3 dispersion correction and all side chains were included. Co−C cleavage is a closedshell to open-shell reaction where bias in the description of unpaired electrons is particularly challenging compared to many inorganic reactions where both sides of the reaction have unpaired electrons. These two situations thus reflect needs for less and more HF exchange, respectively, with 10% as seen in TPSSh representing the “average” case.57,63,64 Thus, the HF exchange “spectrum” should probably be expanded to use 0− 5% for closed-shell to open-shell transitions, and this suggests that the proper choice of HF exchange can be modeled from the change in d-electron occupation of the reaction, as suggested recently by Friesner and co-workers.101 The present results are also consistent with the experience on energy separations between spin states in transition-metal systems57 that also depend critically on the proper balance open- and closed-shell configurations. Energy separations between high-spin and low-spin depend almost linearly on the amount of HF exchange, as first noted by Reiher102 and Paulsen et al.103 and later seen in other cases.89,104,105 This observation gave rise to the construction of B3LYP* with 15% HF exchange.106,107 Thus, a consensus for inorganic chemistry seems to emerge that 10−15% HF exchange provides balanced descriptions of correlation effects in transition-metal systems using current functionals, largely due to its error cancellation of static correlation and self-interaction energies.57

Co−C separation of 3.7 Å in MCAM. This suggests that such a 3.5 Å separation in an enzymatic cage is enough to produce the fully localized Ado radical that then reacts further, lending support to the energetic analysis above that the Co−C separation observed in the enzymes is in fact enough to explain both Co−C labilization and Ado radical formation. The geometry of this state has a fully sp2-hybridized carbon on Ado, as shown in Figure 3C.

■

DISCUSSION The Co−C bond is unique to cobalamins and correspondingly serves a central role in the controlled, reversible production of organic radicals that initiate isomerization reactions on enzymatic substrates.3,5,6 Understanding the nature of this mechanism is of substantial fundamental interest due to its special features and the involvement of a true organometallic alkyl cobalt bond, but also to enable design of catalysts that carry out similar reactions. Thus, the breaking of this bond has been the focus of a large range of experimental and theoretical studies, as mentioned in the introduction of this paper. The present work has investigated the non-enzymatic Co−C cleavage using the largest AdoCbl model so far with all relevant side-chain interactions included quantum-mechanically for the first time. By performing >500 energy calculations of up to 4538 basis functions each and correcting for relativistic, solvent, dispersion, and thermal effects, it yields new, detailed insight into the Co−C cleavage process. With all relevant side chains included, the importance of solvent and dispersion is much larger than previously described, and the accuracy claimed vs experimental data in previous work using small truncated corrins is then almost entirely fortuitous. More precisely, the side chains have substantial dispersion− solvent interactions with the Ado group that affect Co−C cleavage. Still, the full model’s results are not necessarily correct, as multiple approximations prevail, in particular affecting the solvation−dispersion interactions: as the Ado group dissociates, the loss of hydrogen bonding and dispersion interactions with the corrin side chains will be compensated in a real condense phase environment by explicit solvent dispersion terms. The lower side of the Ado group, which packs with the corrin ring (Figure 1) will be solvated, as will the upper side of the corrin that was inaccessible to solvent when Ado was bound. These compensatory interactions are currently impossible to model by quantum chemistry cluster models. This work thus shows that any conclusion based on cluster modeling of a bond dissociation process (or bond-forming process) may be wrong either because the model is unrealistically small or because, as it gets more realistic, dispersion−solvent interactions are poorly described. This is a dilemma of cluster modeling that may be solved by, e.g., reparametrizing solvent−dispersion interactions in the limit of larger, realistic cluster models. However, the detailed quantum mechanics of the Co−C cleavage process is unaffected by solvent dispersion and is a main achievement of the present paper. Notably, the structure of the AdoCbl-dependent enzyme GLUMUT contains two conformations of Ado at ∼3.2 and ∼4.5 Å distance to Co; i.e., the Co−C bond is broken in both conformations.36,100 For MCAM, the Co−C bond is also broken, with Co−C distance of 3.66 Å,90 consistent with these enzymes reducing the equilibrium constant of Co−C cleavage to close to unity. The present computations show that for the full AdoCbl, the Ado radical is fully formed at 3.5 Å, as judged both from sp2-

■

CONCLUSIONS This paper has reported the first quantum-chemical calculations of Co−C cleavage in AdoCbl with full account of the side-chain interactions with the dissociating Ado group, representing the most complete quantum-mechanical treatment of this process so far. The conclusions from the present work are as follow: (i) When one accounts for all interactions within AdoCbl relevant to Co−C dissociation, one finds a very large, combined solvent and dispersion effect acting on this process already in the nonenzymatic reaction. (ii) BP86-D3 remains an accurate method for studying this problem, and B3LYP-D3 displays large underbinding, confirming recently detailed benchmarks for the 7113



■

methylcorrin system82 and showing that claims about the accuracy of B3LYP-D20,21 are erroneous due to comparison of small corrins without side chains to experimental data where the crucial side chains are present. (iii) Ado radical formation, the main deliverable of the AdoCbl system in enzymes, is found to be quantum-mechanically complete at ∼3.5 Å of Co−C dissociation; this explains why crystal structures do not show Co−C bonds shorter than this distance. (iv) BSSEs also contribute to errors (9−12 kJ/mol) in the large realistic AdoCbl model and further favor non-hybrids for these cases. (v) The present work shows that previous studies useing truncated models have been accurate only in cases where the neglect of important side-chain interactions happened to cancel with neglect of dispersion and/or solvent interactions. This is not surprising, given that these side chains directly interact with the large dissociating Ado group. Large quantum-chemical cluster models suffer from challenges relating to the solvent dispersion and sampling problems.18,57 In a condensed phase, the dispersion interactions will be compensated upon dissociation of the large Ado group. This means that the large models, while showing that the small models are incorrect because of loss of specific interactions, are also wrong because the interplay of dispersion, solvent, and BSSE grows with system size. These challenges need to be addressed in future computational chemistry by, e.g., reparametrization of solvent and dispersion models toward larger systems in condensed phases. When expanding the cluster model to the size it has in this work, there is nothing that indicates that the under-binding tendency of B3LYP-D3 will be remedied; in contrast, BP86-D3 still performs most accurately. This is largely attributed to the fact that the Co−C cleavage process is special in having a closed-shell reactant and two open-shell products, which is often not the case in transition-metal chemistry. Such cases seem to require less HF exchange because the open-shell bias only affects the product states. In contrast, reactions staying on one spin surface or even with open-shell singlets will require moreperhaps 10%HF exchange.57

■

REFERENCES

(1) Abeles, R. H.; Dolphin, D. The Vitamin B12 Coenzyme. Acc. Chem. Res. 1976, 9, 114−120. (2) Glusker, J. P. Vitamin B12 and its Coenzymes. Vitamins Hormones 1995, 50, 1−76. (3) Banerjee, R. Radical Carbon Skeleton Rearrangements: Catalysis by Coenzyme B12-Dependent Mutases. Chem. Rev. 2003, 103, 2083− 2094. (4) Banerjee, R.; Ragsdale, S. W. The Many Faces of Vitamin B12: Catalysis by Cobalamin-Dependent Enzymes. Annu. Rev. Biochem. 2003, 72, 209−247. (5) Marsh, E. N. G; Meléndez, G. D. R. Adenosylcobalamin Enzymes: Theory and Experiment Begin to Converge. Biochim. Biophys. Acta 2012, 1824, 1154−1164. (6) Jensen, K. P.; Ryde, U. Cobalamins Uncovered by Modern Electronic Structure Calculations. Coord. Chem. Rev. 2009, 253, 769− 778. (7) Andruniow, T.; Zgierski, M. Z.; Kozlowski, P. M. New Light on the Co−C Bond Activation in B12-Dependent Enzymes from Density Functional Theory. Chem. Phys. Lett. 2000, 331, 509−512. (8) Kozlowski, P. M. Quantum Chemical Modeling of Co−C Bond Activation in B12-Dependent Enzymes. Curr. Opin. Chem. Biol. 2001, 5, 736−743. (9) Freindorf, M.; Kozlowski, P. M. A Combined Density Functional Theory and Molecular Mechanics Study of the Relationship Between the Structure of Coenzyme B12 and its Binding to Methylmalonyl-CoA Mutase. J. Am. Chem. Soc. 2004, 126, 1928−1929. (10) Kuta, J.; Patchkovskii, S.; Zgierski, M. Z.; Kozlowski, P. M. Performance of DFT in Modeling Electronic and Structural Properties of Cobalamins. J. Comput. Chem. 2006, 27, 1429−1437. (11) Dölker, N.; Maseras, F.; Lledos, A. Density Functional Study on the Effect of the trans Axial Ligand of B12 Cofactors on the Heterolytic Cleavage of the Co−C Bond. J. Phys. Chem. B 2003, 107, 306−315. (12) Dölker, N.; Maseras, F.; Lledos, A. A Density Functional Study on the Effect of the Trans Axial Ligand of Cobalamin on the Homolytic Cleavage of the Co-C bond. J. Phys. Chem. B 2001, 105, 7564−7571. (13) Khoroshun, D. V.; Warncke, K.; Ke, S.-C.; Musaev, D. G.; Morokuma, K. Internal Degrees of Freedom, Structural Motifs and Conformational Energetics of 5′-Deoxyadenosyl Radical: Implications for Function in Adenosylcobalamin-Dependent Enzymes. A Computational Study. J. Am. Chem. Soc. 2003, 125, 570−579. (14) Li, X.; Chung, L. W.; Paneth, P.; Morokuma, K. DFT and ONIOM(DFT:MM) Studies on Co−C Bond Cleavage and Hydrogen Transfer in B12-Dependent Methylmalonyl-CoA Mutase. Stepwise or Concerted Mechanism? J. Am. Chem. Soc. 2009, 131, 5115−5125. (15) Jensen, K. P.; Ryde, U. Theoretical Prediction of the Co−C Bond Strength in Cobalamins. J. Phys. Chem. B 2003, 107, 7539−7545. (16) Jensen, K. P.; Ryde, U. Comparison of the Chemical Properties of Iron and Cobalt Porphyrins and Corrins. ChemBioChem 2003, 4, 413−424. (17) Jensen, K. P.; Ryde, U. How the Co−C bond is Cleaved in Coenzyme B12 Enzymes: A Theoretical Study. J. Am. Chem. Soc. 2005, 127, 9117−9128. (18) Ryde, U.; Mata, R. A.; Grimme, S. Does DFT-D Estimate Accurate Energies for the Binding of Ligands to Metal Complexes? Dalton Trans. 2011, 40, 11176−11183. (19) Dölker, N.; Maseras, F.; Siegbahn, P. E. M. Stabilization of the Adenosyl Radical in Coenzyme B12 −A Theoretical Study. Chem. Phys. Lett. 2004, 386, 174−178. (20) Siegbahn, P. E. M.; Blomberg, M. R. A.; Chen, S.-L. Significant van der Waals Effects in Transition Metal Complexes. J. Chem. Theory Comput. 2010, 6, 2040−2044. (21) Chen, S. L.; Blomberg, M. R.; Siegbahn, P. E. M. How is a CoMethyl Intermediate Formed in the Reaction of Cobalamin-Dependent Methionine Synthase? Theoretical Evidence for a Two-Step Methyl Cation Transfer Mechanism. J. Phys. Chem. B 2011, 115, 4066−4077.

ASSOCIATED CONTENT

S Supporting Information *

Electronic energies of the PES in hartrees and kJ/mol for reproduction (Tables S1 and S2); BP86-optimzied PES (Figure S1); details of the effect of Cosmo probe radius (Table S3); vibrational and thermal corrections (Tables S4−S6); energies in hartrees for the BDE calculations (Table S7−S9); relativistic corrections in hartrees (Table S10); structural alignment of Co(II) and Co(III) states (Figure S2); xyz coordinates of optimized structures (Tables S11−S16); ethylene glycol vs water effect on BDE (Table S17); and BDEs computed with various dieclectric constants (Table S18). This material is available free of charge via the Internet at http://pubs.acs.org.

■

Article

AUTHOR INFORMATION

Corresponding Author

*Phone: +45 45 25 24 09. E-mail: [email protected]. Notes

The authors declare no competing financial interest.

■

ACKNOWLEDGMENTS This research has been supported by the Danish Center for Scientific Computing (grant no. 2012-02-23). 7114



Article

(22) Sharma, P. K.; Chu, Z. T.; Olsson, M. H. M.; Warshel, A. A New Paradigm for Electrostatic Catalysis of Radical Reactions in Vitamin B12 Enzymes. Proc. Natl. Acad. Sci. U.S.A. 2007, 104, 9661−9666. (23) Mancia, F.; Keep, N. H.; Nakagawa, A.; Leadlay, P. F.; McSweeney, S.; Rasmussen, B.; Bösecke, P.; Diat, O.; Evans, P. R. How Coenzyme B12 Radicals are Generated: The Crystal Structure of Methylmalonyl-Coenzyme A Mutase at 2 Å Resolution. Structure 1996, 4, 339−350. (24) Marsh, E. N. G. Coenzyme-B12-Dependent Glutamate Mutase. Bioorg. Chem. 2000, 28, 176−189. (25) Babior, B. M.; Carty, T. J.; Abeles, R. H. The Mechanism of Action of Ethanolamine Ammonia-Lyase, A B12-Dependent Enzyme. The Reversible Formation of 5′-Deoxyadenosine from Adenosylcobalamin during the Catalytic Process. J. Biol. Chem. 1974, 249, 1689− 1695. (26) Toraya, T. Radical Catalysis of B12 Enzymes: Structure, Mechanism, Inactivation, and Reactivation of Diol and Glycerol Dehydratases. Cell. Mol. Life Sci. 2000, 57, 106−127. (27) Johnson, J. E., Jr.; Reyes, F. E.; Polaski, J. T.; Batey, R. T. B12 Cofactors Directly Stabilize an mRNA Regulatory Switch. Nature 2012, 492, 133−138. (28) Marsh, E. N. G.; Drennan, C. L. Adenosylcobalamin-Dependent Isomerases: New Insights into Structure and Mechanism. Curr. Opin. Chem. Biol. 2001, 5, 499−505. (29) Hay, B. P.; Finke, R. G. Thermolysis of the Co-C Bond of Adenosylcobalamin. 2. Products, Kinetics, and Co−C Bond Dissociation Energy in Aqueous Solution. J. Am. Chem. Soc. 1986, 108, 4820−4829. (30) Finke, R. G.; Hay, B. P. Thermolysis of Adenosylcobalamin: A Product, Kinetic, and Cobalt−Carbon (C5′) Bond Dissociation Energy Study. Inorg. Chem. 1984, 23, 3041−3043. (31) Brown, K. L.; Zou, X. Thermolysis of Coenzymes B12 at Physiological Temperatures: Activation Parameters for Cobalt-Carbon Bond Homolysis and a Quantitative Analysis of the Perturbation of the Homolysis Equilibrium by the Ribonucleoside Triphosphate Reductase from Lactobacillus leichmannii. J. Inorg. Biochem. 1999, 77, 185− 195. (32) Koenig, T.; Finke, R. G. The Cage Effect and Apparent Activation Parameters for Bond Homolysis. J. Am. Chem. Soc. 1988, 110, 2657−2658. (33) Koenig, T. W.; Hay, B. P.; Finke, R. G. Cage Effects in Organotransition Metal Chemistry: Their Importance in the Kinetic Estimation of Bond Dissociation Energies in Solution. Polyhedron 1988, 7, 1499−1516. (34) Licht, S. S.; Lawrence, C. C.; Stubbe, J. Thermodynamic and Kinetic Studies on Carbon−Cobalt Bond Homolysis by Ribonucleoside Triphosphate Reductase: The Importance of Entropy in Catalysis. Biochemistry 1999, 38, 1234−1242. (35) Reitzer, R.; Gruber, K.; Jogl, G.; Wagner, U. G.; Bothe, H.; Buckel, W.; Kratky, C. Glutamate Mutase from Clostridium cochlearium: The Structure of a Coenzyme B12-Dependent Enzyme Provides New Mechanistic Insights. Structure 1999, 7, 891−902. (36) Gruber, K.; Reitzer, R.; Kratky, C. Radical Shuttling in a Protein: Ribose Pseudorotation Controls Alkyl-Radical Transfer in the Coenzyme B12 Dependent Enzyme Glutamate Mutase. Angew. Chem., Int. Ed. 2001, 40, 3377−3380. (37) Hay, B. P.; Finke, R. G. Thermolysis of the Co−C Bond in Adenosylcorrins. 3. Quantification of the Axial Base Effect in Adenosylcobalamin by the Synthesis and Thermolysis of Axial BaseFree Adenosylcobinamide. Insights into the Energetics of EnzymeAssisted Cobalt−Carbon Bond Homolysis. J. Am. Chem. Soc. 1987, 109, 8012−8018. (38) Marsh, E. N. G.; Holloway, D. E.; Chen, H.-P. In Vitamin B12 and the B12 Proteins; Kräutler, B., Arigoni, D., Golding, B. T., Eds.; Wiley-VCH: Weinheim, Germany, 1998. (39) Vlasie, M. D.; Banerjee, R. Tyrosine 89 Accelerates Co−Carbon Bond Homolysis in Methylmalonyl-CoA Mutase. J. Am. Chem. Soc. 2003, 125, 5431−5435.

(40) Kohn, W.; Becke, A. D.; Parr, R. G. Density Functional Theory of Electronic Structure. J. Phys. Chem. 1996, 100, 12974−12980. (41) Parr, R. G.; Yang, W. Density Functional Theory of Atoms and Molecules; Oxford University Press: New York, 1989. (42) Andruniow, T.; Zgierski, M. Z.; Kozlowski, P. M. Density Functional Theory Analysis of Stereoelectronic Properties of Cobalamins. J. Phys. Chem. B 2000, 104, 10921−10927. (43) Andruniow, T.; Zgierski, M. Z.; Kozlowski, P. M. Theoretical Determination of the Co−C Bond Energy Dissociation in Cobalamins. J. Am. Chem. Soc. 2001, 123, 2679−2680. (44) Jensen, K. P.; Sauer, S. P. A.; Liljefors, T.; Norrby, P.-O. Theoretical Investigation of Steric and Electronic Effects in Coenzyme B12 Models. Organometallics 2001, 20, 550−556. (45) Navizet, I.; Perry, C. B.; Govender, P. P.; Marques, H. M. Cis Influence in Models of Cobalt Corrins by DFT and TD-DFT Studies. J. Phys. Chem. B 2012, 116, 8836−8845. (46) Dybala-Defratyka, A.; Paneth, P. Theoretical Evaluation of the Hydrogen Kinetic Isotope Effect on the First Step of the Methylmalonyl-Coa Mutase Reaction. J. Inorg. Biochem. 2001, 86, 681−689. (47) Jaworska, M.; Lodowski, P. Electronic Spectrum of Co-Corrin Calculated with the TDDFT Method. J. Mol. Struct.: THEOCHEM 2003, 631, 209−223. (48) Kurmaev, E. Z.; Moewes, A.; Ouyang, L.; Randaccio, L.; Rulis, P.; Ching, W. Y.; Bach, M.; Neumann, M. The Electronic Structure and Chemical Bonding of Vitamin B12. Europhys. Lett. 2003, 62, 582− 587. (49) Stich, T. A.; Buan, N. R.; Brunold, T. C. Spectroscopic and Computational Studies of Co2+Corrinoids: Spectral and Electronic Properties of the Biologically Relevant Base-On and Base-Off Forms of Co2+Cobalamin. J. Am. Chem. Soc. 2004, 126, 9735−9749. (50) Liptak, M. D.; Brunold, T. C. Spectroscopic and Computational Studies of Co1+Cobalamin: Spectral and Electronic Properties of the “Superreduced” B12 Cofactor. J. Am. Chem. Soc. 2006, 128, 9144− 9156. (51) Kozlowski, P. M.; Kuta, J.; Galezowski, W. Reductive Cleavage Mechanism of Methylcobalamin: Elementary Steps of Co−C Bond Breaking. J. Phys. Chem. B 2007, 111, 7638−7645. (52) Jensen, K. P.; Ryde, U. Conversion of Homocysteine to Methionine by Methionine Synthase: A Density Functional Study. J. Am. Chem. Soc. 2003, 125, 13970−13971. (53) Jensen, K. P.; Ryde, U. Chemical Properties of Tetrapyrrole Systems with Fe, Co, and Ni. J. Porphyrins Phthalocyanines 2005, 9, 581−606. (54) Kobylianskii, I. J.; Widner, F. J.; Kräutler, B.; Chen, P. Co−C Bond Energies in Adenosylcobinamide and Methylcobinamide in the Gas Phase and in Silico. J. Am. Chem. Soc. 2013, 135, 13648−13651. (55) Becke, A. D. Density-Functional Thermochemistry. III. The Role of Exact Exchange. J. Chem. Phys. 1993, 98, 5648−5652. (56) Lee, C.; Yang, W.; Parr, R. G. Development of the Colle-Salvetti Correlation-Energy Formula into a Functional of the Electron Density. Phys. Rev. B 1988, 37, 785−789. (57) Kepp, K. P. Consistent Descriptions of Metal−Ligand Bonds and Spin-Crossover in Inorganic Chemistry. Coord. Chem. Rev. 2013, 257, 196−209. (58) Frenking, G.; Frohlich, N. The Nature of the Bonding in Transition-Metal Compounds. Chem. Rev. 2000, 100, 717−774. (59) Harvey, J. N. On the Accuracy of Density Functional Theory in Transition Metal Chemistry. Annu. Rep. Prog. Chem. Sect. C 2006, 102, 203−226. (60) Becke, A. D. Density-Functional Exchange-Energy Approximation with Correct Asymptotic Behavior. Phys. Rev. A 1988, 38, 3098−3100. (61) Perdew, J. P. Density-Functional Approximation for the Correlation Energy of the Inhomogeneous Electron Gas. Phys. Rev. B 1986, 33, 8822−8824. (62) Jensen, K. P.; Roos, B. O.; Ryde, U. Performance of Density Functionals for First Row Transition Metal Systems. J. Chem. Phys. 2007, 126, No. 014103. 7115



Article

(63) Jensen, K. P. Bioinorganic Chemistry Modeled with the TPSSh Density Functional. Inorg. Chem. 2008, 47, 10357−10365. (64) Jensen, K. P. Metal− Ligand Bonds of Second-and Third-Row d-Block Metals Characterized by Density Functional Theory. J. Phys. Chem. A 2009, 113, 10133−10141. (65) Tao, J.; Perdew, J. P.; Staroverov, V. N.; Scuseria, G. E. Climbing the Density Functional Ladder: Nonempirical Meta−Generalized Gradient Approximation Designed for Molecules and Solids. Phys. Rev. Lett. 2003, 91, No. 146401. (66) Perdew, J. P.; Tao, J.; Staroverov, V. N.; Scuseria, G. E. MetaGeneralized Gradient Approximation: Explanation of a Realistic Nonempirical Density Functional. J. Chem. Phys. 2004, 120, 6898− 6911. (67) Kepp, K. P. The Ground States of Iron(III) Porphines: Role of Entropy−Enthalpy Compensation, Fermi Correlation, Dispersion, and Zero-Point Energies. J. Inorg. Biochem. 2011, 105, 1286−1292. (68) Kepp, K. P. O2 Binding to Heme is Strongly Facilitated by NearDegeneracy of Electronic States. ChemPhysChem 2013, 14, 3551− 3558. (69) Kepp, K. P.; Dasmeh, P. Effect of Distal Interactions on O2 Binding to Heme. J. Phys. Chem. B 2013, 117, 3755−3770. (70) Ahlrichs, R.; Bär, M.; Häser, M.; Horn, H.; Kölmel, C. Electronic Structure Calculations on Workstation Computers: The Program System Turbomole. Chem. Phys. Lett. 1989, 162, 165−169. (71) Bouquiere, J. P.; Finney, J. L.; Lehman, M. S.; Lindley, P. F.; Savage, H. F. J. High-Resolution Neutron Study of Vitamin B12 Coenzyme at 15 K: Structure Analysis and Comparison with the Structure at 279 K. Acta Crystallogr. 1993, B49, 79−89. (72) Ouyang, L.; Rulis, P.; Ching, W. Y.; Nardin, G.; Randaccio, L. Accurate Redetermination of the X-ray Structure and Electronic Bonding in Adenosylcobalamin. Inorg. Chem. 2004, 43, 1235−1241. (73) Weigend, F.; Ahlrichs, R. Balanced Basis Sets of Split Valence, Triple Zeta Valence and Quadruple Zeta Valence Quality for H to Rn: Design and Assessment of Accuracy. Phys. Chem. Chem. Phys. 2005, 7, 3297−3305. (74) Klamt, A.; Jonas, V.; Bürger, T.; Lohrenz, J. C. W. Refinement and Parametrization of COSMO-RS. J. Phys. Chem. A 1998, 102, 5074−5085. (75) Jensen, K. P. Computational Studies of Modified [Fe3S4] Clusters: Why Iron is Optimal. J. Inorg. Biochem. 2008, 102, 87−100. (76) Andrade, S. G.; Gonçalves, L. C. S.; Jorge, F. E. Scaling Factors for Fundamental Vibrational Frequencies and Zero-Point Energies Obtained from HF, MP2, and DFT/DZP and TZP Harmonic Frequencies. J. Mol. Struct.: THEOCHEM 2008, 864, 20−25. (77) Grimme, S.; Antony, J.; Ehrlich, S.; Krieg, H. A Consistent and Accurate Ab Initio Parametrization of Density Functional Dispersion Correction (DFT-D) for the 94 Elements H-Pu. J. Chem. Phys. 2010, 132, 154104. (78) Blomberg, M. R. A.; Borowski, T.; Himo, F.; Liao, R.-Z.; Siegbahn, P. E. M. Quantum Chemical Studies of Mechanisms for Metalloenzymes. Chem. Rev. 2014, 114, 3601−3658. (79) Hirao, H. Which DFT Functional Performs Well in the Calculation of Methylcobalamin? Comparison of the B3LYP and BP86 Functionals and Evaluation of the Impact of Empirical Dispersion Correction. J. Phys. Chem. A 2011, 115, 9308−9313. (80) Grimme, S. Accurate Description of van der Waals Complexes by Density Functional Theory Including Empirical Corrections. J. Comput. Chem. 2004, 25, 1463−1473. (81) Grimme, S. Semiempirical GGA-Type Density Functional Constructed with a Long-Range Dispersion Correction. J. Comput. Chem. 2006, 27, 1787−1799. (82) Kozlowski, P. M.; Kumar, M.; Piecuch, P.; Li, W.; Bauman, N. P.; Hansen, J. A.; Lodowski, P.; Jaworska, M. The Cobalt−Methyl Bond Dissociation in Methylcobalamin: New Benchmark Analysis Based on Density Functional Theory and Completely Renormalized Coupled-Cluster Calculations. J. Chem. Theory Comput. 2012, 8, 1870−1894. (83) Perdew, J. P.; Burke, K.; Ernzerhof, M. Generalized Gradient Approximation Made Simple. Phys. Rev. Lett. 1996, 77, 3865−3868.

(84) Schafer, A.; Horn, H.; Ahlrichs, R. Fully Optimized Contracted Gaussian Basis Sets for Atoms Li to Kr. J. Chem. Phys. 1992, 97, 2571− 2577. (85) Kozlowski, P. M.; Kamachi, T.; Kumar, M.; Yoshizawa, K. Reductive Elimination Pathway for Homocysteine to Methionine Conversion in Cobalamin-Dependent Methionine Synthase. J. Biol. Inorg. Chem. 2012, 17, 611−619. (86) Rovira, C.; Kozlowski, P. M. First Principles Study of Coenzyme B12. Crystal Packing Forces Effect on Axial Bond Lengths. J. Phys. Chem. B 2007, 111, 3251−3257. (87) Furche, F.; Perdew, J. P. The Performance of Semilocal and Hybrid Density Functionals in 3d Transition-Metal Chemistry. J. Chem. Phys. 2006, 124, 044103. (88) Zhao, Y.; Truhlar, D. G. Comparative Assessment of Density Functional Methods for 3d Transition-Metal Chemistry. J. Chem. Phys. 2006, 124, 224105. (89) Jensen, K. P.; Cirera, J. Accurate Computed Enthalpies of Spin Crossover in Iron and Cobalt Complexes. J. Phys. Chem. A 2009, 113, 10033−10039. (90) Froese, D. S.; Kochan, G.; Muniz, J. R.; Wu, X.; Gileadi, C.; Ugochukwu, E.; Krysztofinska, E.; Gravel, R. A.; Oppermann, Y. U.; Yue, W. W. Structures of the Human GTPase MMAA and Vitamin B12-Dependent Methylmalonyl-CoA Mutase and Insight into Their Complex Formation. J. Biol. Chem. 2010, 285, 38204−38213. (91) Cowan, R. D.; Griffin, D. C. Approximate Relativistic Corrections to Atomic Radial Wave Functions. J. Opt. Soc. Am. 1976, 6, 1010−1014. (92) Jensen, K. P.; Ryde, U. The Axial N-base has Minor Influence on Co−C Bond Cleavage in Cobalamins. J. Mol. Struct.: THEOCHEM 2002, 585, 239−255. (93) Kräutler, B.; Keller, W.; Kratky, C. Coenzyme B12 Chemistry: The Crystal and Molecular Structure of Cob(II)alamin. J. Am. Chem. Soc. 1989, 111, 8936−8938. (94) Xhang, Y.; Cremer, P. S. Interactions between Macromolecules and Ions: the Hofmeister Series. Curr. Opin. Chem. Biol. 2006, 10, 658−663. (95) Jansen, H. B.; Ros, P. Non-Empirical Molecular Orbital Calculations on the Protonation of Carbon Monoxide. Chem. Phys. Lett. 1969, 3, 140−143. (96) Van Duijneveldt, F. B.; Van Duijneveldt-Van de Rijdt, J. G. C. M.; Van Lenthe, J. H. State of the Art in Counterpoise Theory. Chem. Rev. 1994, 94, 1873−1885. (97) Boys, S. F.; Bernardi, F. The Calculation of Small Molecular Interactions by the Differences of Separate Total Energies. Some Procedures with Reduced Errors. Mol. Phys. 1971, 19, 553−566. (98) Pyykko, P. Relativistic Effects in Structural Chemistry. Chem. Rev. 1988, 88, 563−594. (99) Chowdhury, S.; Banerjee, R. Thermodynamic and Kinetic Characterization of Co−C Bond Homolysis Catalyzed by Coenzyme B12-Dependent Methylmalonyl-CoA Mutase. Biochemistry 2000, 39, 7998−8006. (100) Gruber, K.; Kratky, C. Coenzyme B12 Dependent Glutamate Mutase. Curr. Opin. Struct. Biol. 2002, 6, 598−603. (101) Hughes, T. F.; Friesner, R. A. Correcting Systematic Errors in DFT Spin-Splitting Energetics for Transition Metal Complexes. J. Chem. Theory Comput. 2011, 7, 19−32. (102) Reiher, M. Theoretical Study of the Fe(phen)2(NCS)2 SpinCrossover Complex with Reparametrized Density Functionals. Inorg. Chem. 2002, 41, 6928−6935. (103) Paulsen, H.; Duelund, L.; Winkler, H.; Toftlund, H.; Trautwein, A. X. Free Energy of Spin-Crossover Complexes Calculated with Density Functional Methods. Inorg. Chem. 2001, 40, 2201−2203. (104) Zein, S.; Borshch, S. A.; Fleurat-Lessard, P.; Casida, M. E.; Chermette, H. Assessment of the Exchange-Correlation Functionals for the Physical Description of Spin Transition Phenomena by Density Functional Theory Methods: All the Same? J. Chem. Phys. 2007, 126, No. 014105. (105) Daku, L. M. L.; Vargas, A.; Hauser, A.; Fouqueau, A.; Casida, M. E. Assessment of Density Functionals for the High-Spin/Low-Spin 7116



Article

Energy Difference in the Low-Spin Iron(II) Tris(2,2′-bipyridine) Complex. ChemPhysChem 2005, 6, 1393−1410. (106) Reiher, M.; Salomon, O.; Hess, B. A. Reparameterization of Hybrid Functionals based on Energy Differences of States of Different Multiplicity. Theor. Chem. Acc. 2001, 107, 48−55. (107) Paulsen, H.; Trautwein, A. X. Density Functional Theory Calculations for Spin Crossover Complexes. Top. Curr. Chem. 2004, 235, 197−219.

7117


CoâC Dissociation of Adenosylcobalamin (Coenzyme B12): Role of

Recommend Documents

CoâC Dissociation of Adenosylcobalamin (Coenzyme B12): Role of