Unexpected Routes of the Mutagenic Tautomerization of the T Nucleobase in the Classical A·T DNA Base Pairs: A QM/QTAIM Comprehensive View

In this paper using quantum-mechanical (QM) calculations in combination with Bader's quantum theory of “Atoms in Molecules” (QTAIM) in the continuum with ε = 1, we have theoretically demonstrated for the first time that revealed recently highly-energetic conformers of the classical A·T DNA base pairs – Watson-Crick [A·T(wWC)], reverse Watson-Crick [A·T(wrWC)], Hoogsteen [A·T(wH)] and reverse Hoogsteen [A·T(wrH)] – act as intermediates of the intrapair mutagenic tautomerization of the T nucleobase owing to the novel tautomerisation pathways: A·T(wWC)↔A·T*(w⊥WC); A·T(wrWC)↔A·TO2*(w⊥rWC); A·T(wH)↔A·T*(w⊥H); A·T(wrH)↔A·TO2*(w⊥rH). All of them occur via the transition states as tight ion pairs (A+, protonated by the N6H2 amino group)·(T−, deprotonated by the N3H group) with quasi-orthogonal geometry, which are stabilized by the participation of the strong (A)N6+H···O4−/O2−(T) and (A)N6+H···N3−(T) H-bonds. Established tautomerizations proceed through a two-step mechanism of the protons moving in the opposite directions along the intermolecular H-bonds. Initially, proton moves from the N3H imino group of T to the N6H2 amino group of A and then subsequently from the protonated N6+H3 amino group of A to the O4/O2 oxygen atom of T, leading to the products – A·T*(w⊥WC), A·TO2*(w⊥rWC), A·T*(w⊥H), and A·TO2*(w⊥rH), which are substantially non-planar, conformationally-labile complexes. These mispairs are stabilized by the participation of the (A)N6H/N6H'···N3(T) and (T)O2H/O4H···N6(A) H-bonds, for which the pyramidalized amino group of A is their donor and acceptor. The Gibbs free energy of activation of these mutagenic tautomerizations lies in the range of 27.8–29.8 kcal·mol−1 at T = 298.15 K in the continuum with ε = 1.


INTRODUCTION
Clarification at the microstructural level of the physico-chemical mechanisms underlying the formation of the mutagenic tautomers of the DNA bases via the mutagenic tautomerization of the classical Watson-Crick DNA base pairs is a matter of extreme importance for such branches of life science as molecular biophysics and molecular biology, since it enables us to understand the sources of the genome instability (Watson and Crick, 1953a,b;Löwdin, 1963Löwdin, , 1966Topal and Fresco, 1976). Genome instability is frequently associated with mutations in DNA, playing role in cancer development due to DNA replication errors (Liu et al., 2014;Tomasetti et al., 2017).
Mutagenic tautomerization of the DNA bases attracts researchers' curiosity since the establishment of the spatial architecture of DNA molecule (Watson and Crick, 1953a) and further formulation of the tautomeric hypothesis of the origin of the spontaneous point mutations by Watson and Crick (Watson and Crick, 1953b).
Distinguished quantum chemist Per-Orlov Löwdin proposed original idea based on the electronic structure of the complementary A·T and G·C pairs of the DNA bases (Löwdin, 1963(Löwdin, , 1966, which makes possible their conversion into the high-energy tautomerized states -A * ·T * (L) and G * ·C * (L) base pairs [currently known as Löwdin's base pairs; here and below rare, in particular mutagenic (Brovarets' and Hovorun, 2010a;, tautomers are marked with an asterisk] causing origin of the transitions and transversions during the DNA replication. Löwdin believed that these transformations should be carried out by the double proton transfer (DPT) in the opposite directions along the neighboring intermolecular hydrogen (H) bonds through the quantum tunneling. These representations played an extremely important role in the formation of new visions in quantum biology and attracted the attention of a wide range of Löwdin's followers (Florian et al., 1994;Gorb et al., 2004;Bertran et al., 2006;Cerón-Carrasco and Jacquemin, 2013;Maximoff et al., 2017).
However, from the physico-chemical point of view it was established that generally accepted Löwdin's mechanism of the DPT along the intermolecular H-bonds in the Watson-Crick DNA base pairs cannot be the source of the formation of the mutagenic tautomers of the nucleobases due to the absence of the reverse barrier of tautomerization in the A·T(WC) pair of the DNA bases and its small value in comparison with kT (0.62 kcal·mol −1 at T = 298.15 K) for the G·C(WC) DNA base pair (Gorb et al., 2004;Bertran et al., 2006;Brovarets' et al., 2012;Hovorun, 2014a,b, 2015a).
Recently, we have proposed another mechanism of the mutagenic tautomerization of the A·T(WC) and G·C(WC) pairs of the DNA bases, which is alternative to Löwdin's approach, occurring via the sequential intrapair proton transfer and shifting of the bases relative each other, which ultimately leads to the wobble configuration (Brovarets' and Hovorun, 2015b). Moreover, we have discovered this intrinsic ability to perform wobble↔Watson-Crick / Watson-Crick↔wobble tautomeric transitions via the sequential intrapair proton transfer for all possible incorrect base mispairs, which are active players in the field of the spontaneous point mutagenesis: purine·pyrimidine -G·T and A·C (Brovarets' and Hovorun, 2009, 2015c, purine·purine -A·A, A·G and G·G (Brovarets' and Hovorun, 2015e,f) and pyrimidine·pyrimidine -C·C, C·T and T·T (Brovarets' and Hovorun, 2015f,g). Notably, these interconverisons are accompanied by a significant rebuilding of the base mispairs with Watson-Crick architecture into the mismatches wobbled toward both minor and major DNA grooves and vice versa. Moreover, it was established that these tautomerisation reactions occur non-dissociatively and are accompanied by the consequent replacement of the unique patterns of the intermolecular specific interactions along intrinsic reaction coordinate (IRC) (Brovarets' et al., 2013(Brovarets' et al., , 2017a.
These data allows to suggest that the intrapair tautomeric transition of the wobble pairs from the main tautomeric form into the rare, mutagenic, having a WC or close to its configuration, and vice versa, is the key to understanding of the microstructural mechanisms of the emergence of the spontaneous transitions and transversions at the DNA replication (Brovarets' and Hovorun, 2009, 2015b,c,d, 2016. Moreover, these theoretical approaches have been partly experimentally confirmed for some DNA/RNA purine·pyrimidine pairs (Nedderman et al., 1991(Nedderman et al., , 1993Kimsey et al., 2015Kimsey et al., , 2018. In this study, we succeeded to further elaborate such approach and to reveal new mechanism of the mutagenic tautomerization of the classical A·T DNA base pairs (Scheme 1) as their intrinsic property, lying beyond classical representations at the microstructural level and which was not presented in the literature before. For the first time, it was theoretically shown using QM/QTAIM methods, that the transition of these pairs into the substantially non-planar, high-energy conformers  provokes intrapair mutagenic tautomerization of the T DNA base from the canonical, diketo into the rare, enol tautomeric forms T * and T * O2 Hovorun, 2014a, 2015b,d;. Moreover, for the first time we have investigated in details conformationally-tautomeric properties of the classical A·T DNA base pairs (Brovarets' et al., 2018b,c,d,e).
Transition states (TSs) of these mutagenic tautomerisations are tight ion pairs (A + , protonated by the N6H 2 amino group; T − , deprotonated by the N3H group) with quasi-orthogonal geometry, which are stabilized by the participation of the strong (A)N6 + H· · · O4 − /O2 − (T) and (A)N6 + H· · · N3 − (T) H-bonds. Discovered reaction of the mutagenic tautomerization proceeds through the stepwise mechanism of the PT along the H-bonds: primarily proton moves from the imino group N3H of T to the N6H 2 amino group of A and then proton transfers from the protonated N6 + H 3 amino group of A to the O4/O2 oxygen atom of T, leading to the products, which are substantially nonplanar, conformationally-labile complexes. These complexes are stabilized by the participation of the (A)N6H/N6H ′ · · · N3(T) and (T)O2H/O4H· · · N6(A) H-bonds, for which the pyramidalized amino group of A DNA base acts as their donor and acceptor. The Gibbs free energy of the activation of the mutagenic tautomerizations lies in the range of 27.79-29.83 kcal·mol −1 at T = 298.15 K in the continuum with ε = 1.

COMPUTATIONAL METHODS
Geometries of the investigated DNA base pairs and TSs of their mutual tautomeric and conformational transformations, as well as their harmonic vibrational frequencies were calculated at the B3LYP/6-311++G(d,p) level of theory (Hariharan and Pople, 1973;Krishnan et al., 1980;Lee et al., 1988;Parr and Yang, 1989;Tirado-Rives and Jorgensen, 2008), using Gaussian'09 package (Frisch et al., 2009) followed by the IRC calculations in the forward and reverse directions from each TS using Hessian-based predictor-corrector integration algorithm (Hratchian and Schlegel, 2005). A scaling factor that is equal to 0.9668 Hovorun, 2010b,c,d, 2011;El-Sayed et al., 2015) was applied in this study for the correction of the harmonic frequencies of all DNA base pairs and TSs of their tautomeric and conformational transitions. We have confirmed the TSs, localized by Synchronous Transit-guided Quasi-Newton method (Peng et al., 1996), on the potential energy landscape by the presence of one and only one imaginary frequency in the vibrational spectra of the complexes. We applied standard TS theory for the estimation of the activation barriers of the tautomeric transformations (Atkins, 1998). Single point electronic energy calculations have been performed using MP2 level of theory (Frisch et al., 1990) and aug-cc-pVDZ Dunning's cc-type basis set (Kendall et al., 1992), which was confirmed as appropriate level of theory for the analogous systems and tasks (Lozynski et al., 1998;Danilov et al., 2005;Matta, 2010;Rutledge and Wetmore, 2012;Pérez-Sánchez, 2016, 2017;. All calculations were performed for the base pairs in the continuum with a dielectric constant of ε = 1 as their intrinsic property, that is adequate for modeling of the processes occurring in real systems (Bayley, 1951;Dewar and Storch, 1985;Petrushka et al., 1986;García-Moreno et al., 1997;Mertz and Krishtalik, 2000;Bebenek et al., 2011;Wang et al., 2011;Maximoff et al., 2017) without deprivation of the structurally functional properties of the bases in the composition of DNA Pérez-Sánchez, 2016, 2017;.
The Gibbs free energy G for all structures was obtained in the following way: where E el -electronic energy, while E corr -thermal correction. The Gibbs free energy of activation or barrier for the forward tautomeric/conformational transition was calculated as the difference between the Gibbs free energy of the TS and reactant of the reaction. The Gibbs free energy for the reverse tautomeric/conformational transition was calculated as the difference between the Gibbs free energy of the TS and product of the reaction.
Electronic interaction energies E int were calculated at the MP2/6-311++G(2df,pd) level of theory as the difference between the total energy of the base pair and energies of the monomers and corrected for the basis set superposition error (BSSE) (Boys and Bernardi, 1970;Gutowski et al., 1986) through the counterpoise procedure (Sordo et al., 1988;Sordo, 2001).
Bader's quantum theory of Atoms in Molecules (QTAIM) (Bader, 1990;Matta and Hernández-Trujillo, 2003;Matta, 2014;Lecomte et al., 2015) was applied to analyse the electron density distribution, using software package AIMAll (Keith, 2010). The presence of the bond critical point (BCP), namely (3,−1) BCP, and a bond path between hydrogen donor and acceptor or between two electronegative covalently bonded atoms, as well as the positive value of the Laplacian at this BCP ( ρ > 0), were considered as criteria for the H-bond or attractive van der Waals contact formation (Matta et al., 2006;Hovorun, 2014c, 2018b;. Wave functions were obtained at the level of theory used for geometry optimisation. The energies of the attractive van der Waals contacts (Matta and Boyd, 2007; in the TSs of the conformational transitions of the tautomerized base pairs were calculated by the empirical Espinosa-Molins-Lecomte (EML) formula (Espinosa et al., 1998;Mata et al., 2011), based on the electron density distribution at the (3,−1) BCPs of the specific contacts: in this formula V(r) is a value of a local potential energy at the (3,−1) BCP. The energies of the conventional AH···B H-bonds were evaluated by the empirical Iogansen's formula (Iogansen, 1999): in this formula ν is a magnitude of the frequency shift of the stretching mode of the AH H-bonded group involved in the AH···B H-bond relatively the unbound group. The partial deuteration was applied in order to avoid the effect of vibrational resonances (Brovarets' and Hovorun, 2015h;. The atomic numbering scheme for the DNA bases was conventional (Saenger, 1984).

OBTAINED RESULTS AND DISCUSSION
In our previous study, for the first time we have succeeded to establish in the classical biologically-important A·T Frontiers in Chemistry | www.frontiersin.org DNA base pairs with C s symmetry -Watson-Crick (WC), reverse Watson-Crick A·T(rWC), Hoogsteen A·T(H) and reverse Hoogsteen A·T(rH) DNA base pairs (Scheme 1) (Donohue and Trueblood, 1960;Haschemeyer and Sobell, 1963;Hoogsteen, 1963;Brovarets', 2013a,b;Yang et al., 2015;Poltev et al., 2016;Zhou, 2016;Szabat and Kierzek, 2017) -novel high-energetic, dynamically-stable, mirrorsymmetrical A·T(w WC ) R,L , A·T(w H ) R,L , A·T(w rWC ) R,L and A·T(w rH ) R,L conformational states (Figure 1) . Their distinguished feature is significantly non-planar structure (C 1 symmetry), which is caused by the pyramidal structure of the ≥C6N6H 2 amino fragment of the A DNA base, which amino group acts simultaneously as a donor and an acceptor of the specific intermolecular interactions with T DNA base by two (T)N3H· · · N6(A) and (A)N6H/N6H ′ · · · O4/O2(T) Hbonds (the N6H ′ bond has trans-orientation relatively the N1C6 bond of A). Each of the four A·T Watson-Crick DNA base pairs transfers into the aforementioned conformers via two mirror-symmetric pathways through the TS A·T(WC)↔A·T(w WC )R,L , TS A·T(rWC)↔A·T(w rWC )R,L , TS A·T(H)↔A·T(w H )R,L and TS A·T(rH)↔A·T(w rH )R,L (C 1 symmetry). At this, mirror-symmetrical complexes, which are enantiomers, are marked with the subscripts R and L. Notably, enantiomers in the achiral environment demonstrate identical scalar physicochemical characteristics and differ only by the direction of the dipole moment.
Possible biological role of these conformers was also elucidated, in particular -their participation in the nondissociative conformational interconversions of all four classical A·T DNA base pairs (Brovarets' et al., 2018b,e). Recently, we have identified novel pathway of the mutagenic tautomerisation of these structures through the quasi-orthogonal transition state as A − · T + (Brovarets' et al., 2018c).
These data inspired us to elaborate further this novel point of view for the classical objects such as biologically-important A·T DNA base pairs and allow to suggest the possibility of the mutagenic tautomerization of T through the stepwise PT along the appropriate intermolecular H-bonds from the N3H imino group of T to the N6 atom of the N6H 2 amino group of A in the just-mentioned conformers and further -from the protonated amino group NH + 3 of A to the O4/O2 oxygen atoms of T depending on the starting pair.
It was established that novel pathways of the mutagenic tautomerization of the T DNA base in the classical A·T DNA base pairs (Scheme 1) are initiated by their spontaneous conformational transition into the high-energy A·T(w WC ) R,L , A·T(w H ) R,L , A·T(w rWC ) R,L and A·T(w rH ) R,L conformers as well as are controlled by the TSs as tight ion pairs (A + , protonated by the N6H 2 amino group)·(T − , deprotonated by the N3H imino group) with electronic energy of interaction E int ∼145 (19.00) and TS A + ·T − A·T(w rH )R,L↔A·T * O2 (w ⊥ rH )L,R (21.48 kcal·mol −1 ) are characterized by the quasi-orthogonal arrangement of the bases relatively each other and are stabilized by the participation of the two non-equivalent strong H-bonds (A)N6 + H· · · N3 − (T) and (A)N6 + H· · · O4 − /O2 − (T) [the first of them is significantly weaker (∼15.1-18.6 kcal·mol −1 ), than the second one (∼22.5-23.1 kcal·mol −1 )]. Protonated amino group N6 + H 3 of A for these TSs acts simultaneously as donor and acceptor of the Hbonding and has such spatial orientation, that its N6 + H/N6 + H ′ bond, which is not involved in the H-bonding with T, lies in the plane of the purine ring (Figure 1, Tables 1, 2).
It is worth to mention that each of the investigated tautomeric and conformational transitions proceed through two mirrorsymmetric pathways and do not change cys/trans mutual orientation of the N1H and N9H glycosydic bonds of the bases. At the mutagenic tautomeric transformations of the DNA bases some R/L structures transfer into the other L/R structures and vice versa (Figures 1, 2).
Terminal tautomerized complexes are conformationally-labile and pairwise interconvert into each other according to four mechanisms (Tables 1, 2).
Two of these tautomerization reactions are controlled by the TSs -TS 1,2 A·T * (w ⊥ WC )R,L↔A·T * (w ⊥ H )R,L (14.0, 10.9 cm −1 ) and TS 1,2 A·T * O2 (w ⊥ rWC )R,L↔A·T * O2 (w ⊥ rH )R,L (13.4, 11.1 cm −1 ) with low values of the imaginary frequencies provided in the brackets. At this, one-single intermolecular (T)O4H/O2H· · · N6(A) Hbond between the O4H/O2H hydroxyl groups of T * /T * O2 and N6 nitrogen atom of the piramidalized amino group of A participates in the stabilization of the TS 1 s. In the case of TS 2 s, when T hangs over A, the (T)O4H/O2H· · · N6(A) H-bond coexists together with attractive van der Waals contacts with significantly increased ellipticity -N3· · · C6 and O2· · · C4 in the case of TS 2 A·T * (w ⊥ WC )↔A·T * (w ⊥ H ) and N3· · · C6 in the case of TS 2 A·T * O2 (w ⊥ WC )↔A·T * O2 (w ⊥ H ) ( Table 2). Notably, conformational transformations, which are controlled by the TS 1 s are the most energetically favorable (1.86 and 1.92) in comparison with the TS 2 s (2.56 and 2.63 kcal·mol −1 ) ( Table 1). In these cases R/L structures are converted into the other R/L structures.
All tautomeric and conformational transitions without exceptions are dipole-active processes, since they are TABLE 1 | Energetic characteristics (in kcal·mol −1 ) of the discovered mutagenic tautomerizations of the T DNA base in the classical A·T DNA base pairs via the DPT and conformational transformations of their products obtained at the MP2/aug-cc-pVDZ//B3LYP/6-311++G(d,p) level of QM theory in the continuum with ε = 1 at T = 298.15 K (see Figures 1, 2).   (Iogansen, 1999), EML (Espinosa et al., 1998;Mata et al., 2011; marked with an asterisk) or Nikolaienko-Bulavin-Hovorun (Nikolaienko et al., 2012; marked with a double asterisk) formulas, kcal·mol −1 . h The dipole moment of the complex, D.
Frontiers in Chemistry | www.frontiersin.org accompanied by a noticeable change in the dipole moment of the involved complexes (Table 2).
Interestingly, that among all without exception investigated in this work H-bonded structures, the total energy of the intermolecular specific contacts (H-bonds and attractive van der Waals contacts) contribute only a part of the electron energy of the monomer interactions (0.26-0.98; see Figures 1, 2). This result is in a good agreement with the previously published data for the others H-bonded pairs of nucleotide bases (Brovarets' and Hovorun, 2014d).
Notably, the methyl group of the T DNA base does not change its orientation during all, without exception, processes of the tautomeric and conformational transformations. Moreover, the heterocycles of the DNA bases remain planar, despite their ability for the out-of-plane bending (Govorun et al., 1992;Hovorun et al., 1999;Nikolaienko et al., 2011).
Finally, we would like to emphasize the fact that the presence of the conformational transitions between the complexes -products of the A·T * (w ⊥ WC )R, L↔ A·T * (w ⊥ H )R, L and A·T * O2 (w ⊥ rWC )R, L↔ A·T * O2 (w ⊥ rH )R, L, tautomerizations indicating the close structural relationship between tautomerization the classical A·T(WC) and A·T(H) DNA base pairs, on the one hand, and A·T(rWC) and A·T(rH), on the other hand (Brovarets' et al., 2018b,e).

CONCLUSIONS
In this study, we came out from the existing framework of the mechanisms of the origin of the mutagenic tautomerization of the classical A·T DNA base pairs (Brovarets', 2013b;Brovarets' et al., 2018a,b,c,d,e).
Here we have shed light on the revealed for the first time physico-chemical mechanism of the intrapair mutagenic tautomerization of the T DNA base within the novel highly-energetic conformers of the classical A·T DNA base pairs -Watson-Crick [A·T(w WC )], reverse Watson-Crick [A·T(w rWC )], Hoogsteen [A·T(w H )] and reverse Hoogsteen [A·T(w rH )], which have been analyzed in details in our previous paper . These reactions -A·T(w WC )↔A·T * (w ⊥ WC ), A·T(w rWC )↔A·T * O2 (w ⊥ rWC ), A·T(w H )↔A·T * (w ⊥ H ), A·T(w rH )↔A·T * O2 (w ⊥ rH ) proceed through the stepwise proton transfer via the TSs as tight A + ·T − ion pairs, which Gibbs free energy of activation lies in the range of 27.79-29.83 kcal·mol −1 at T=298.15 K, thus creating the substantially nonplanar, conformationally-labile complexes -A·T * (w ⊥ WC ), A·T * O2 (w ⊥ rWC ), A·T * (w ⊥ H ) and A·T * O2 (w ⊥ rH ). Furthermore, formed complexes involving mutagenic T * /T * O2 tautomers are able to conformationally interconvert between each other according to reaction pathways -A·T * (w ⊥ WC )↔A·T * (w ⊥ H ) and A·T * O2 (w ⊥ rWC )↔A·T * O2 (w ⊥ rH ).

AUTHOR CONTRIBUTIONS
OB, study conception and design, acquisition of data, drafting of manuscript analysis and interpretation of data, performance of calculations, discussion of the obtained data, preparation of the numerical data for Tables, graphical materials