1. Introduction

The study of enzyme mechanisms is of vast importance and has become central to both bioorganic chemistry and computational chemistry. Studies in this field by Bruice, Benkovic, Jencks, and Bender over the past 40 years have contributed largely to understanding the mode and scope by which certain enzymes achieve their catalytic activities [1]. These and studies of many others have been carried out in order to understand enzymatic catalysis which is characterized by high substrate specificity, chemoselectivity, and stereospecificity along with large rate enhancements. Among these are (1) the proximity model of Bruice [2, 3] on the intramolecular cyclization of dicarboxylic semiesters, (2) the “orbital steering” theory suggested by Koshland [4, 5], (3) the “spatiotemporal hypothesis” devised by Menger [6–9] which describes the importance of the distance between the two reactive centers in determining whether a reaction is inter- or intramolecular, and (4) the gem-trimethyl lock (stereopopulation control) proposal of Cohen [10–12].

Reaction models for mimicing enzyme catalysis usually fail to reach a desirable rate enhancement due to a lack of a high degree of conformational restrictions comparable to that existing in the enzyme-substrate complex [1]. In 1970, Cohen studied the lactonization of a series of hydroxyhydrocinnamic acids and found accelerations above 10¹⁵ when compared to the intermolecular cyclization of the corresponding counterparts. He attributed this large enhancement to what he called “stereopopulation control” [10–12]. Cohen’s interpretation of his data was criticized by various researchers who claimed that the high rate enhancement results from relief in strain energy during the lactonization of the hydroxyhydrocinnamic acid and not because of a stereopopulation control within the trimethyl lock system [13–15].

Recently, we have been engaged in exploring the driving force(s) behind the remarkable acceleration rates in some intramolecular reactions [16–19]. Using ab initio molecular orbital at different levels, molecular mechanics, and semiempirical molecular orbital methods, we studied the thermodynamic and kinetic behavior of the lactonization of some hydroxy acids, the cyclization reactions of Bruice’s di-carboxylic semiesters, the proton transfer reaction in Menger’s system, and the cyclization of some -bromoalkanecarboxylate anions. The results from these studies revealed the following main conclusions. (1) Rate enhancement in intramolecular reactions can be driven by proximity orientation which is due to strain effects or by proximity that is not related to strain effects of a starting material and/or a corresponding transition state. For example, our study on cyclization of Bruice’s dicarboxylic semiesters reveals that the activation energy in these systems is dependent on the difference in the strain energies of the transition states and the reactants, and there is no relationship between the cyclization rate and the distance between the nucleophile and the electrophile, and the reactivity extent of the semiester system is linearly correlated with the strain energy difference between the transition state and the reactant. On the other hand, acid lactonization of hydroxy acids reveals that the enhancement in rates of the lactonization in this system stems from the close proximity of the electrophile to the nucleophile. Further, it shows that the rate of the lactonization reaction is solely dependent on the ratio between the angle of attack of the nucleophile and the distance between the two reacting centers. This is in accordance with Menger’s “spatiotemporal hypothesis” that relates distance between the nucleophile and the electrophile to the rate of the reaction. (2) Significant rate enhancements in intramolecular reactions are due to both enthalpic and entropic effects and not only due to enthalpic effects as was proposed by Bruice. (3) The nature of the reaction (intermolecular or intramolecular) is largely dependent on the distance between the two reacting centers. For example, our ab initio calculations on Menger’s system show that when the distance between the two reacting centers is 2.4 Å, the reaction is intramolecular, whereas, when the distance is 3 Å, the reaction prefers the intermolecular process. Further, our study shows that the proximity between the nucleophile and electrophile is largely dependent on the strain energy of the system. For a strained system, the distance between the two reacting centers is shorter than that in unstrained systems [16–19].

In this letter, we describe the DFT at B3LYP/6-31G () and the RHF at 6-31G levels, as well as AM1 semiempirical calculations results (thermodynamics and kinetics) for the acid-catalyzed and uncatalyzed lactonization reactions of the trimethyl lock system. Our goal was to establish the driving force behind the remarkable accelerations of the intramolecular reaction in the tri-methyl lock system 1c in Figure 1.

Figure 1: Hydroxyhydrocnnamc acids.

2. Methods

The DFT, HF, and AM1 calculations were carried out using the quantum chemical package Gaussian-98 [20]. The MM2 molecular mechanics strain energy calculations were performed using Allinger’s MM2 program installed in Chem 3D Ultra 8.0 [21]. The starting geometries of all the molecules in this study were obtained using the Argus Lab program [22]. The calculations were carried out based on the restricted Hartree-Fock (RHF) method with full optimization of all geometrical variables (bond lengths, bond angles, and dihedral angles) [23]. To avoid results with local minima optimization, the keyword Freq Opt = (Z-matrix, MaxCycle = 300, CalcAll) GFINPUT IOP(6/7 = 3) was used in the input files of the starting geometries. The geometry optimizations included estimations of second derivatives (Hessian matrix) for each of the parameters in each species ( for planar structures) [24]. DEP analytical gradients were used throughout the optimization. Geometries were optimized in internal coordinates and were terminated when Herbert’s test was satisfied in the Broyden-Fletcher-Goldfarb-Shanno (BFGS) method.

An energy minimum (a stable compound or a reactive intermediate) has no negative vibrational force constant. A transition state is a saddle point which has only one negative vibrational force constant [25].The transition state structures were verified by their only one negative frequency. The verification was accomplished by viewing the frequency results via the Molden program [26]. The “reaction coordinate method” [27–29] was used to calculate the activation energy for the lactonization processes of hydroxyl acids 1a–1d. In this method, the value of one bond is limited for the appropriate degree of freedom while all other variables are optimized. The activation energy values for the cyclization reactions were calculated from the difference in the energies of the global optimum structures for the reactants 1a–1d and the derived transition states of the cyclization reactions. The transition state structures for the cyclization reactions of 1a–1d were obtained from the decrease in the distance between the phenolic oxygen (O10) and the carbonyl carbon (C1) in increments of 0.1 Å. Full optimization of the transition states was accomplished after removing the constrains imposed while executing the energy profile. The frequency results obtained from the optimization were viewed by Molden program and it was found that all the transition state structures, studied here, have only one negative frequency. The DFT at B3LYP/ 6-31G () and HF at 6-31G levels of the reactions of 1a–1d were calculated with and without the inclusion of solvent (water, dielectric constant = 78.39). The keywords SCF = Tight and SCRF = Dipole were used in the input files when calculating energies with the incorporation of a solvent.

3. Results and Discussion

In order to determine whether stereopopulation control, suggested by Cohen et al. [10–12] or conventional relief of strain energy, proposed by Winans and Wilcox [13] and supported by a theoretical study by Houk et al. [15] is the driving force for the enzyme-like accerelation in 1c, we have calculated the AM1 thermodynamic properties for the lactonization of hydroxyhydrocinnamic acids 1a–1d (Figure 1) using Gaussian 98 version 3.0 [20] available at our Al-Quds computer center. These thermodynamic calculations were performed (the calculated data is depicted in Tables S1 and S2, Supplementary Material, available at doi:10.1155/2009/240253) in order to prove that rates of reactions cannot be predicted from free energy changes (MM2 strain energy, thermodynamic parameter) as was previously reported by Winans and Wilcox [13].

The AM1 calculations of the enthalpic and entropic energies of the lactonization processes of hydroxy acids having a trimethyl lock (1c) and lacking the trimethyl lock (1a) (see Table 1) reveal that the difference in the free energy values between the lactonization reactions of the two systems () is 9.7 kcal/mol ( = tetrahedral intermediate- hydroxyl acid). Based on this value, the calculated rate enhancement of the lactonization of 1c compared to that of 1a should be 1.1 10⁷, which is 10⁴-fold less than the experimentally determined value (see Table 1) [10–12]. Further, the MM2 calculated difference in strain energies (, for the numbering see Figure 1) in the two systems 1a and 1c is about 7 kcal/mol, which is equal to 1.2 10⁵ (10⁶-fold less than the experimentally value). This result excludes the notion that the remarkable rate enhancement in the lactonization of 1c is solely accommodated by conventional relief of strain as was reported by the groups of Houk and Wilcox [13–15].

Table 1: AM1 and MM2 calculated thermodynamic properties of hydroxyhydrocinnamic acid derivatives. (All values were calculated by AM1 method except for the values of Es which were calculated by MM2 method.)

To test whether the conformational restriction plays an important role in the rate accerelation of the lactonization processes of hydroxy acids 1a and 1c, calculations of the rotational barrier around the C2–C3 bond were executed. Note that the C1–O10 distance in hydroxyhydrocinnamic acids is dependent upon the C1/C2/C3/C4 dihedral angle. This C1–O10 distance should be crucial for enhancing the lactonization rate according to the stereopopulation control (conformational locking) suggested by Cohen [10–12]. When plotting the C1–O10 distance in 1a and 1c against the heats of formation of the resulting conformations (see Figure 2), an important result emerges. In 1c, the shortest C1–O10 distance values correspond to the most stable conformations, whereas for 1a the shortest distance values correspond to the highest energy conformations.

Figure 2: (2a) a plot of energy versus dihedral angle C1C2C3C4 in 2a and 2c. (2b) a plot of energy versus C1–O10 distance in 2a and 2c.

Further, calculations of the thermodynamic parameters (, and ) of each of the energetically most stable conformers in 1a, as well as those for the conformer with the shortest C1–O10 distance, reveal that for 1a to reside in the most productive position to react intramolecularly requires a = 10.73 kcal/mol (composed of = 4.91 kcal/mol and - = 5.82 kcal/mol). In contrast, 1c requires = 3.93 kcal/mol ( = 1.85 kcal/mol and- = 2.08 kcal/mol) in order to fulfill the same task (see Table S3, Supplementary Material, available at doi:10.1155/2009/240253). The roughly 7 kcal/mol difference in ground-state free energy for the two systems is equivalent to about 1.4 10⁶ in rate enhancement of 1c over 1a. Again, the calculations show clearly that the conformational locking is not the critical factor in the rate accerelation seen with 1c (comprising only a calculated 10⁶ value versus an experimental 10¹¹ value) [10–12].

In light of the results of the AM1 calculations that show neither the steric effects suggested by Winans and Wilcox [13] nor the conformational locking suggested by Cohen [10–12] are the main driving force for the enhancement in the lactonization of 1c, we have calculated, using the “reaction coordinate” method, the activation energy values (^‡) for the formation and the collapse of the tetrahedral intermediate involved in the lactonization process. The DFT at B3LYP/6-31G () and HF/6-31G levels as well as the AM1 activation energy values were calculated with and without the inclusion of solvent (water, dielectric constant 78.39) and the results obtained indicate that the effect of water on the relative rate values is negligible (1.2–1.8 kcal/mol for the three different lactonization processes of 1a, 1c, and 1d). This is in accordance with previously reported studies of Houk et al. on lactonization of hydroxy acids, that indicate that the solvation effect more-or-less cancels out when comparing reactivities of species having the same structural features (even though the absolute rate constants cannot be evaluated) [13–15]. Thus, reaction rates were computed from the activation energy results for both the acid-catalyzed and uncatalyzed lactonization processes of hydroxy acids 1a–1d. The calculated rates are linearly related to the experimentally determined rate values, the correlations which are shown in (1) (for the acid-catalyzed process) and in (2) (for the un-catalyzed process). The AM1 calculations results of the activation energy values for the lactonization processes of hydroxy acids 1a–1d were correlated very well with the values obtained using the ab initio method at the DFT B3LYP/6-31G () and RHF/6-31G, and the correlation results are depicted in (3)-(4), respectively: where is the cyclization rate of hydroxy acid 1b-1d, and is the cyclization rate of hydroxy acid 1a; where is the cyclization rate of protonated hydroxy acid 2b-2d, and is the cyclization rate of protonated hydroxyl acid 2a; where (HF) and (AM1) are the calculated activation energy values by RHF/6-31G and AM1 methods, respectively; where (B3LYP) is the calculated activation energy values by B3LYP/6-31G ().

The calculated rate values of the lactonization processes of 1a and 1c using (2) are 10⁴ and 10¹⁴, respectively. This gives a predicted rate enhancement of about 10¹⁰ which is in good agreement with the experimentally determined value (10¹¹). Our calculations also reveal that the rate-limiting step is formation rather than collapse of the tetrahedral intermediate [16–19], in opposition to the conclusions reported by Houk et al. [15].

It is accepted that strain-accelerated reactions occur for compounds that, by necessity, are rigid with high bond-rotation barriers that exceed those required for the reaction. Further, such compounds have distorted bond distances or/and bond angles when compared to strain-free compounds [1]. Our theoretical calculations indeed predict distortion of bond angles and bond distances in the phenyl ring of 1c which is supported by the X-ray crystal structure of the corresponding alcohol [30]. However, the same distortions with the same magnitude are observed with 1b and to some extent with the corresponding tetrahedral intermediates 3c and 3b. If the acceleration is due to strain relief, we should see comparable rates for 1c and 1b, but actually the lactonization rate of 1c is 10¹⁰ times faster than that of 1b. This suggests no significant strain relief upon the lactonization of 1c.

The second convincing point excluding strain as the main driving force for the rate acceleration is that in 1c the rotation barrier around C2-C3 is found to be 3 kcal/ mol smaller than that of 1a. This surprising conclusion arises from the fact that 1c has a stronger intramolecular carboxyl/hydroxyl hydrogen-bonding than does 1a. If Winans and Wilcox were correct [13] and the rate acceleration is indeed due to strain-relief, then the C2-C3 bond in 1a should be more affordable to rotation than that of 1c. In other words, if strain-accelerated reactions occur in compounds that are rigid with high bond rotation barriers, then the rate comparison between 1a and 1c cannot be due to strain effects.

For a better understanding the lactonization process, we conducted reaction-coordinate calculations for 1a and 1c (when the C1–O10 distance is in the range of  Å) by calculating the change in energy as a function of the O10/C1/C2 attack angle () and the C1–C10 distance between the two reacting centers (). An excellent correlation was observed between the entalpic energy and the ratio between the attack angle and the distance. Further, it was found that the slope follows the order . This result suggests that the energy needed to increase the value of angle α to reach the optimal value for formation of the transition state is less in the case of 1c when compared to 1a. This in turn indicates that the approach of the hydroxyl group to the carbonyl carbon in the case of 1c is much easier than in the case of 1a.

4. Conclusions

The combined results indicate that the driving force in the first and second stages of the approach (C1–O10 =  Å, and  Å), is due to a proximity effect. Hydroxy acids that are rigid in the organized state have lower activation energies (such as 1c) than those with less rigidity (such as 1a). It is worthy to note that Mengerhas advocated abandoning thermodynamic models involving entropy in favor of distance effects on rate, and in fact he has derived an equation relating rate and distance [31, 32].

Acknowledgments

The author thanks the Karamans Co. and the German-Palestinian-Israeli fund agency for support of our hardware computational facilities. Special thanks go to Dr. Omar Deeb and Sherin Alfalah and Donia Karaman for support in computational software and technical assistance. His sincere appreciations are given to Professor Fred Menger (Emory University, Ga, USA) for helpful discussions.