Thus, considering the above data on the influence of substituents, polarization of the substrate and hydrogen bonding we cannot clearly identify the key factor for facilitating the reaction. Neither the substrate structure itself nor the expected H-bonding interactions can describe the enhancement of reaction rate in the [4+2] cycloaddition reaction.

The calculated influence of hydrogen bonding revealed a possible change in the activation barrier in the range of 2–5 kcal/mol in a good agreement with previous studies [ 28 , 37 ]. Such a change could contribute to the reaction rate enhancement, although it can not be suggested as the only leading factor for accelerating the cycloaddition reaction under experimental conditions [ 14 , 38 ].

As shown by the density functional calculations, a coordination of Gln to the hydroxyl group has a small influence on the reaction barrier (<1 kcal/mol); on the PM6 energy surface the reaction became slightly more favored by ~3 kcal/mol (entries 1, 2; Table 2 ). The influence of the second amino acid slightly lowered the barrier at the B3LYP and O3LYP levels, while the opposite change was found for the M062X, PM6 and PM6-DH2 energy surfaces (entries 3, 2; Table 2 ). A coordination of Ser to the carbonyl group of the substrate slightly lowered the activation energy (entries 4, 3; Table 2 ), which was in line with proposed polarization in the substrate ( Fig. 2 ). At the density functional level the activation barriers were decreased by 0.4–5.7 kcal/mol and by 1.3 kcal/mol at the PM6 level ( Table 2 ). A replacement of hydrogen bonding of amino acids by explicit solvation with water molecules resulted in a similar change in the calculated parameters (ΔE ≠ 2-TS = 33.8 kcal/mol at the PM6 level).

The third important point is a possible influence of hydrogen bonding on the reactivity of the substrate in the cycloaddition reaction. These interactions may facilitate cycloaddition in the active site of enzymes due to modulation of energies of the HOMO/LUMO molecular orbitals by specific hydrogen bonding. Positioning of a hydrogen bond acceptor to diene (increasing the energy of HOMO) and a hydrogen bond donor to alkene (decreasing the energy of LUMO) could narrow the energy gap and lower the activation barrier. Both functional groups of amino acids and water molecules from a solvent shell may be involved in hydrogen bonding with the substrate. It was demonstrated that a coordination of glutamine or asparagine amino acids with diene and serine or tyrosine with dienophile can lower the activation barrier by 4.7 kcal/mol [ 37 ]. We verified the role of both factors and studied the influence of coordination of one, two and three amino acids, as well as interactions with water molecules. The calculations were carried out at various density functional levels (B3LYP, O3LYP, M062X) and semiempirical level (PM6) including an accurate description of H-bonding and dispersion interactions (PM6-DH2). PM6-DH2 has shown the same accuracy as sophisticated B2-PLYP-D/TZVPP and M06-2X/6-311+G(3df,2p) levels. The deviation of 0.4 kcal/mol with benchmark CCSD(T)/CBS level was re-ported [ 47 , 48 ].

The second factor of interest is the degree of polarization of the dienophylic part of the substrate and charge delocalization over the conjugated chain ( 2-TS and 2-TS’ ; Fig. 2 ) proposed as a plausible channel to overcome the transition state [ 14 , 38 ]. An analysis of the charge delocalization showed insignificant alterations upon the formation of the transition state (see S2 Table ; S3 Table ; S4 Table ). Thus, the calculations did not confirm that the zwitterionic (or significantly polarized) nature of the transition state would necessarily facilitate the reaction.

The calculated activation barriers and reaction energies at the B3LYP level showed only minor changes in spite of significant structural variations ( Table 1 ). A rather small difference of 1.5 kcal/mol was found for reactions a and e . The stepwise assembly of the target substrate maintained the activation barriers in a narrow range of ΔE ≠ = 23.6–26.1 kcal/mol for the [4+2] cycloaddition reactions they are involved in. The calculated ΔH energy surface showed the same trend ( Table 1 ). Thus, the decomposition analysis demonstrated that the considered structural changes neither led to a dramatic decrease of the activation barrier in the [4+2] cycloaddition in the biosynthesis of Spinosyn A, nor could account for carrying out the reaction under mild room temperature conditions. The calculated ΔG free energy surface showed an expected decrease in the activation barrier and reaction energy due to an appearance of the first (cf. a and b ) and second (cf. d and e ) intramolecular linkages. This effect is well-known and reflects lowering of the reaction entropy [ 39 , 40 ]. In overall, the calculated energy parameters agree well with other [4+2] cycloaddition reactions involving 1,3-diene and alkene units as well as with an experimental study of reaction a [ 42 , 43 , 44 , 45 , 46 ].

First, we carried out a chemical structure decomposition analysis by identifying principal structural fragments and revealing their influence on energetic parameters. The following key-stages were involved in the decomposition analysis ( Fig. 2 ): a) the reference [4+2] cycloaddition reaction; b) attachment of an intramolecular-(CH 2 ) 3 - linkage; c) attachment of a hydroxyl group in the linkage; d) elongation of the conjugated chain; e) attachment of a carbonyl group and macrocyclic linkage. In fact, the structure decomposition follows an increasing complexity of assembly of Spinosyn A over the Diels-Alder reaction core in a step-by-step manner. Full geometry optimization of initial structures, transition states and products were performed for each of the reactions ( Fig. 2 a— 2e ). The concerted transition states were successfully located for all the cases studied (see S1 Fig ; S2 Fig ; S3 Fig ; S4 Fig ; S1 Table and S1 File for structures and geometry parameters). The increase of molecular complexity led to a more pronounced asynchrony of the C-C bonds formation. The result is in a good agreement with the earlier DFT calculations of the [4+2] transition states [ 41 , 42 , 43 , 44 , 45 , 46 , 83 ].

Structural decomposition analysis of the cycloaddition reaction involved in the biosynthesis of Spinosyn A into principal components (the atom numbering was the same as in compound 1 for comparative purpose).

The classical Diels-Alder reaction involves cycloaddition of butadiene and ethylene units and results in a cyclohexene ring via a single transition state ( Fig. 2a ) [ 1 , 2 , 3 , 4 , 8 , 9 , 10 , 11 , 12 ]. In the biosynthesis of Spinosyn A, several structural changes were introduced to this [4+2] reaction core, and an important question is to what extent these changes may affect the reactivity of the ring construction. Possible contributions of the following known factors have been considered in the present study: 1) influence of substituents in the substrate; 2) polarization of the substrate; and 3) influence of hydrogen bonding.

An alternative route for the studied [4+2] cycloaddition reaction

Upon considering hydrogen bonding structures we have revealed another mechanistic pathway of enzymatic catalysis of the cycloaddition reaction. The course of the Diels-Alder reaction is governed by shortening of intramolecular distances resulting in the formation of C-C bonds, which is reflected by an increase of potential energy of the system until reaching the transition state. In our finding, the energetic profile of the studied system may adopt a more favorable curve if the movement along the reaction coordinate is accompanied by stepwise hydrogen bonding. Particularly, a computational study of the stepwise process involving i) H-bonding followed by shortening of intramolecular distances; ii) the next H-bonding followed by further shortening of intramolecular distances; and iii) repeat the same sequence until complete reaction; revealed a much more favorable reaction profile.

In the studied system this kind of transformation may start with initial state I, where substrate 1 is located in a hydrophobic environment of the enzyme (Fig. 3). For the substrate bond distances C 7 -C 11 = 3.217 Å and C 4 -C 12 = 3.952 Å were calculated at the starting point (see Fig. 2 for atom numbering). A coordination of the substrate with an amino acid residue initiates an arrangement of the enzyme active center and H-bond formation with C 9 -OH was calculated to be exothermic by 7.3 kcal/mol (II; Fig. 3) (we have selected Gln and Ser for the model computational study as these amino acids were shown to be involved in the active center of enzymes with cycloaddition activity (see refs. [35, 36]) and the choice of these amino acids is not strictly limited, since coordination of other amino acids would lead to a similar energy gain). Upon moving in the enzyme active center the substrate experienced the influence of the protein environment and returned to the thermoneutral state on the calculated energy surface (III) after contraction of the macrocycle resulted in a shortening of the C 7 -C 11 and C 4 -C 12 bonds to 2.619 Å and 2.830 Å, respectively.

On the next step, H-bond formation between the C 17 -OH group and the second amino acid led to the energy gain of 12.7 kcal/mol (IV), and made further contraction of the forming bonds to C 7 -C 11 = 2.190 Å and C 4 -C 12 = 2.690 Å possible. At this point the calculated energy of the studied system was again around the thermoneutral point (V).

The third step involved H-bonding of the C 15 = O group (VI) and a macrocycle ring contraction to the bond distances C 7 -C 11 = 2.010 Å and C 4 -C 12 = 2.610 Å. The process furnished structure VII with the relative energy around zero point (Fig. 3); however, geometric parameters of VII were significantly shifted towards the “product-like” state as compared to initial structure I.

Here we described the process involving movement along the reaction coordinate in an H-bonding followed by ring contraction manner (Fig. 3). Exactly the same energy surface was calculated when the opposite sequence of ring contraction followed by H-bonding was studied. Indeed, both “H-bonding/contraction” and “contraction/H-bonding” type movements along the reaction coordinate led to the same stationary points that can be unambiguously defined by the values of the C1–C6 and C4–C5 bonds lengths.

Starting with VII as an initial structure, the cycloaddition reaction took place easily through the transition state VIII-TS (Fig. 4) and involved overcoming a small calculated activation barrier of 6.3 kcal/mol (Fig. 3). As expected, the cycloaddition reaction is exothermic due to formation of C-C bonds and due to release of internal strain of the contracted molecule.

When moving along the reaction path, the stepwise contraction of the substrate was energetically balanced by the H-bond formation. The substrate motion via the energetically balanced reaction path furnished formation of structure VII, which was more “product-like” as compared to I, while the overall energy of the system changed in the range from 0.0 to -12.7 kcal/mol. Considering the overall surface, the energy difference of 19.0 kcal/mol between the highest and lowest points (IV and VIII-TS) would contribute to the rate-determining activation parameters.

For comparison, cycloaddition of the isolated molecule of I requires overcoming the energy barrier as high as 36.6 kcal/mol (XI-TS; Fig. 3). It is important to note that a simple coordination of the same amino acids (without macrocycle contraction) did not lead to a dramatic decrease in the activation energy (Table 2 and discussion above). An important role in such a process plays an influence of the protein environment in order to govern directed macrocycle ring contraction.

Thus, the catalytic effect of the calculated pathway lowered the activation energy by 17.6 kcal/mol: from 36.6 kcal/mol for XI-TS to 18.0 kcal/mol for VIII-TS. The reaction energy remained the same for both pathways—decoordination of amino acids from IX results in X (Fig. 3). To confirm the reliability of the selected computational approach single point calculations were carried out at CCSD(T)/cc-pVDZ level. The CCSD(T) calculations have clearly indicated significant lowering of the activation energy by 13.0 kcal/mol (17.6 kcal/mol at DFT level) and a similar overall reaction energy -10.8 kcal/mol (-11.9 kcal/mol at DFT level).