- Email: [email protected]
Imidazopyridopyrimidine base pairing motifs consisting of four hydrogen bonds: a quantum chemical study
Chemical Physics Letters 392 (2004) 508–513 www.elsevier.com/locate/cplett
Imidazopyridopyrimidine base pairing motifs consisting of four hydrogen bonds: a quantum chemical study Jirı Czernek
*
Institute of Macromolecular Chemistry, Academy of Sciences of the Czech Republic, Heyrovsky Square 2, 162 06 Praha, The Czech Republic Received 28 January 2004; in final form 16 April 2004 Available online 19 June 2004
Abstract Following a recent experimental report on the synthesis and thermal stability of imidazopyridopyrimidine-containing structures with the ability to form four hydrogen bonds, the geometries and relative stabilities of the pairs of model bases have been studied using the HF, MP2 and B3LYP approaches. The 1 H NMR chemical shifts of the base pairs have been estimated by the SOS-DFPTIGLO method and employed in the analysis of the experimental proton spectra. Ó 2004 Elsevier B.V. All rights reserved.
1. Introduction DNA is the structural support of genetic material and therefore, the key factor in many vital cellular processes. While double-stranded right-handed helix is a regular conformation adopted by DNA in cells, there is an increasing number of results pointing to the biological importance of alternative structures such as bulges, hairpins, cruciforms, branched junctions, quadruplexes, and left-handed Z conformations [1]. Hydrogen bonding (H-bonding) [2] is of paramount importance in stabilizing helical (through the Watson–Crick base pairing [3]) and also the other DNA structures, which involve alternative pairing schemes such as Hoogsteen type interactions [4]. Consequently, a considerable attention has been paid to the design, synthesis and structural analysis of short oligodeoxynucleotides (ODNs) featuring different H-bonding motifs [5]. Most recently, two sets of new base pair combinations consisting of four hydrogen bonds have been proposed (see Fig. 1) and incorporated into ODNs [6]. Importantly, the tricyclic skeletons adopted by Minakawa et al., [6] i.e., imidazo[50 ; 40 :4,5]pyrido[2,3-d] pyrimidine nucleosides, possess the ability to form both Watson–Crick and Hoogsteen types of H-bonding and thus are potentially *
Tel.: +420-296809290; Fax: +420-296809410. E-mail address: [email protected] (J. Czernek).
0009-2614/$ - see front matter Ó 2004 Elsevier B.V. All rights reserved. doi:10.1016/j.cplett.2004.06.007
useful in stabilizing both the double-stranded helical and alternative DNA structures. These novel pairing motifs have been characterized using thermal stabilities of several duplexes of ODNs containing different number of the tricyclic nucleosides. It has been shown that in some cases the thermal stability of imidazopyridopyrimide pairs surpasses that of their Watson–Crick counterparts. Presently, high-level ab initio computational techniques are able to reliably account for the intrinsic interactions in nucleic acids, i.e., to quantitatively describe the nature and energetics of various interactions in (sizeable fragments of) nucleic acids without the presence of external effects such as solvent, crystal packing and others (see [7] for a recent review). The knowledge of the intrinsic molecular interactions is indispensable in understanding the role played by respective structural elements, including intercalators [8] and metal cations [9], in DNA structure and function. Thanks to a significant benchmarking effort, useful guidelines for correct and computationally efficient computational studies of nucleic acids have been devised [10]. As these approaches also provide for accurate geometrical characteristics of the studied structures, an advantage can be taken from recent advances in ab initio methods for calculating NMR (nucleic magnetic resonance) parameters [11] to predict and interpret the spectra of nucleosides [12,13] and base pairs [14–16].
J. Czernek / Chemical Physics Letters 392 (2004) 508–513
(a)
H4b
H2
N
N
O
H4a
H2
N
HN
N
H5
N
NH
H9
N
H6
N
H9 N
N N
O
H7a
H7b
diamino (1)
dioxo (2)
H4b
(b)
H2
N
N HN
N
H9
N
O
H4a
H5
H6
H2
N
N
NH
N
H9 N
N O
H7a
N H7b
amino-oxo (3)
oxo-amino (4)
Fig. 1. Base pairing in diamino–dioxo (1 2,a) and amino-oxo–oxoamino (3 4,b) complexes. The numbering of protons for which the chemical shielding has been evaluated is shown.
The main objective of this Letter, is the quantum chemical description of the structure and stability of the diamino–dioxo and amino-oxo–oxo-amino heterodimers of imidazopyridopyrimide base pairs (Fig. 1). The results are compared with the computational data obtained at the same, sufficiently high levels of theory for guanine–cytosine (GC) and adenine–thymine (AT) Watson–Crick base pairs [17], and for selected dimers containing four hydrogen bonds [18], and are discussed in the context of the thermal stabilities found by Minakawa et al. [6]. Moreover, to aid in the interpretation of the 1 H NMR spectra of different pairing motifs, the amino-oxo and oxo-amino homodimers have been studied in addition to the above-mentioned heterodimers, and the chemical shifts of all the base pairs have been computed using a density functional theory (DFT)based method.
2. Calculations The structures of four dimers formed by imidazo[50 ; 40 :4,5]pyrido[2,3-d]pyrimidine 1 (denoted as ‘diamino’ in [6]), imidazo[50 ; 40 :4,5]pyrido[2,3-d]pyrimi-
509
dine-4,7(5H,6H)-dione 2 (‘dioxo’), imidazo[50 ; 40 :4,5]pyrido[2,3-d]pyrimidin-7(6H)-one 3 (‘amino-oxo’), and imidazo[50 ; 40 :4,5]pyrido[2,3-d]pyrimidin-4(5H)-one 4 (‘oxo-amino’) were studied [19], namely, the heterodimers 1 2 and 3 4, and the homodimers 3 3 and 4 4 (see Fig. 1). Thus, the molecular geometries of these dimers and of the monomers of 1–4 were optimized using the Hartree–Fock (HF) theory and the 6-31G** basis set of atomic orbitals, and the density functional theory with B3LYP (Becke’s three parameter exchange [20] and Lee, Yang, Parr correlation [21]) functionals and the 6-311G** basis set. The Cs symmetry was imposed on the monomers and heterodimers, while the C2h symmetry was exploited in the calculations of the homodimers. To verify the character of the optimized structures and to obtain the data for thermodynamic computations (see below), calculations of the harmonic vibrational frequencies were subsequently performed at the same level as the geometrical optimizations. The interaction energies, DE, were calculated using the variational supermolecular approach as the energy difference between the complex and the isolated components. Namely, for the structures optimized at the HF/6-31G** level, the second-order Møller–Plesset (MP2) energies were obtained in the frozen core approximation using the standard 6-31G basis set augmented with one set of diffuse d-polarization functions with the exponent of 0.25, which were placed on the second-row atoms. This methodology will be referred to as MP2/6-31G*(0.25)// HF/6-31G**. To assess the importance of correlation energy for the stabilization of imidazopyridopyrimide pairs, the HF/6-31G*(0.25) interaction energies were evaluated as well. Moreover, for the structures optimized using the B3LYP/6-311G** approach, the DE’s were computed at the same level of theory. All interaction energies were corrected for the basis set superposition error by means of the counterpoise method of Boys and Bernardi [22]. In addition, for the B3LYP/ 6-311G** data the deformation energy, EDEF , was computed by subtracting the energy of the optimized isolated monomers from the energies of the monomers in the dimer geometry [10]. All the above-mentioned computations were performed using the default algorithms and settings of GA U S S I A N 98 [23]. Employing the interaction energies and the unscaled values of harmonic vibrational frequencies, the B3LYP/6-311G** enthalpies and free energies of the complexes at 298 K were estimated in the rigid rotor-harmonic oscillator-ideal gas approximation [24]. The 1 H NMR chemical shielding tensors were obtained using the SOS-DFPT-IGLO [25] method as implemented in the deMon-MASTER-CS code [26,27]. The SOS-DFPT-IGLO strategy combines the sum-overstates density functional (Rayleigh–Schr€ odinger) perturbation theory with the individual gauges for localized (molecular) orbitals [28]. The IGLO-III basis set of
510
J. Czernek / Chemical Physics Letters 392 (2004) 508–513
Kutzelnigg et al. [28] (roughly of quadruple-n quality with two sets of polarization functions), the Perdew– Wang [29,30] exchange-correlation potential, the FINE angular integration grid with 64 radial shells [31], and the approximation Loc. 1 SOS-DFPT [25] were applied. The molecular orbitals were localized by the method of Boys [32].
3. Results and discussion Table 1 summarizes the counterpoise–corrected interaction energies as provided by the well-tested [10] MP2/6-31G*(0.25)//HF/6-31G** strategy for the 1 2, 3 4, 3 3, and 4 4 dimers. DE’s of an important type of structures which have been studied using the abovementioned approach, i.e., the GC and AT Watson– Crick base pairs [17], are shown for comparison purposes. The data obtained for GC and AT have been shown, due to a compensation of errors, to lie close to the results of the CCSD(T)/aug-cc-pVDZ calculations [33]. Thus, the MP2/6-31G*(0.25)//HF/6-31G** approach is expected to yield reliable characteristics also of the imidazopyridopyrimide bases. Data in Table 1 clearly indicate a significant stabilizing interaction in all the dimers investigated, with DE’s ranging from 91.3 (in the case of 3 3) to 136.0 kJ/mol (1 2). These values are comparable to the complexation energy of the GC pair, but they are ca. two times higher than DE of the AT Watson–Crick base pair (see Table 1). It has been demonstrated that the stabilization of hydrogen-bonded bases of nucleic acids [10], and of similar types of hydrogen-bonding systems [34], is predominantly electrostatic in origin. Hence, quite reliable characteristics of this type of structures can be obtained already at the HF level with medium-sized polarized basis sets [33]. This holds also for four imidazopyridopyrimide base pairs studied here: the order of the HF/6-31G*(0.25) interaction energies is identical to that found by the MP2 method. The correlation contribution to the stabiliza-
Table 1 The negative of total MP2 interaction energy, DEMP2 , and its HF, DEHF , and electron correlation, DECOR , components (kJ/mol), calculateda for different complexes Complex
DEMP2
DEHF
DECOR
12 34 33 44 GC AT
136.0 98.0 105.5 91.3 108.4 52.1
117.8 78.7 85.0 72.9 103.3 40.7
18.1 19.3 20.5 18.4 5.0 11.3
See the text for details. a The results for the GC and AT base pairs are taken from Table 1 of [17].
tion energy, i.e., the difference between the respective MP2 and HF data, constitutes 20% of the DEMP2 in the case of 3 4, 3 3 and 4 4 pairs, and 13% in the case of the most stable dimer 1 2. This contribution amounts to 5% and 22% for the GC and AT pairs, respectively (Table 1). It is noted that all vibrational frequencies obtained for the HF/6-31G** geometries have been found to be real except the one reflecting inversion of the amino group [10] in the position 7 of the diamino monomer 1. Let us now turn to the discussion of the B3LYP/ 6-311G** results. This approach, which now can be routinely applied to quite sizeable systems, is known to be highly accurate. For example, it has been shown to be capable of providing smaller absolute deviations for the bond lengths and angles than the QCISD/6-31G* method [35], and successfully applied in a recent comparative study of the geometries and relative stabilities of 17 multiply hydrogen-bonded complexes [18]. Consequently, the B3LYP/6-311G** technique was adopted here to study the geometry, energetics, and thermodynamics of the imidazopyridopyrimide base pairs. As was the case of the HF/6-31G** geometries, all B3LYP/ 6-311G** vibrational frequencies are real except of the one describing inversion of the amino group in the position 7 of the compound 1. The data relevant to the thermodynamics of imidazopyridopyrimide dimers formation are summarized in Table 2 together with the results for the most and the least stable complexes containing four hydrogen bonds, taken from the set recently studied by Lukin and Leszczynski, [18] i.e., the structures denoted as 11 and 14, respectively (see also Fig. 2). In accord with the MP2/6-31G*(0.25)//HF/6-31G**computed interaction energies, a significant stabilization of imidazopyridopyrimide dimers is predicted on the basis of the B3LYP/6-311G** calculations. These two approaches also agree in the ordering of the DE’s for respective pairs, with the 3 3 being the least and the 1 2 the most stable dimers (their DE’s amount to 96.1 and Table 2 The negative of the DFT-computed: interaction energy, DEDFT , deformation energy, EDEF , and standard molar enthalpy, DH 0 298 , entropy, DS 0 298 , and Gibbs energy, DG0 298 , changes at 298 K calculateda for different complexes (kJ/mol)b Complex
DEDFT
EDEF
DH 0 298
DS 0 298
DG0 298
12 34 33 44 11 14
146.7 102.5 110.0 96.1 196.4 58.1
)17.5 )9.9 )11.2 )8.9 NA NA
119.5 82.9 91.6 75.3 191.3 69.7
188.2 184.4 201.6 173.5 249.3 213.8
63.3 27.9 31.5 23.6 117.0 6.0
See the text for details. The results for the complexes denoted as 11 and 14 are taken from Table 2 of [18]. b DS 0 298 is given in J/(mol K). a
J. Czernek / Chemical Physics Letters 392 (2004) 508–513
H O
N O
N
H
H
N
N
N
N
H
H
N
O
CH3
H3C O
N H
11 O H
H N H O
N H N CH3
N
N
H
H
N
O NH
14 Fig. 2. The structure of the complexes denoted as 11 and 14.
146.7 kJ/mol accordingly). In agreement with previous findings [10] is the fact that the larger the stabilization, the higher is the (always repulsive) deformation energy of the bases (see Table 2). An inclusion of the zero-point vibration energy and temperature-dependent enthalpy contributions moderately reduces the stabilization (as expressed by the values of DH 00 and DH 0298 accordingly). The entropic term is fairly uniform (16% variation within the series), but significant, constituting up to 69% of the standard molar enthalpy change at 298 K. However, its inclusion does not alter the stability order of imidazopyridopyrimide base pairs as expressed by their DG0298 ’s (Table 2). Clearly, apart from the solvent effects, there are many structural factors affecting the stability of hydrogenbonded systems, including the type of donor–acceptor sequences, specific substituents, and conformational properties. Nevertheless, it can be of some interest to compare the results obtained for imidazopyridopyrimide base pairs with those calculated at the B3LYP/ 6-311G** level for some complexes containing four hydrogen bonds [18]. The energetic and thermodynamic data describing the most and the least stable assemblies studied by Lukin and Leszczynski [18] generally encompass the values found for imidazopyridopyrimide dimers. However, the exception is the (destabilizing) entropic term, which is for imidazopyridopyrimides bases smaller than for both the dimers 11 and 14 (see
511
Fig. 2 and Table 2). In the case of the least stable structure 14 from [18], the enthalpy stabilization is almost cancelled out by the entropy contribution, leading to the standard Gibbs energy change at 298 K of only 6.0 kJ/mol. On the other hand, the structure 11 is significantly more stable than imidazopyridopyrimide base pairs, partly due to the presence of two intramolecular hydrogens bonds [18]. Minakawa et al. (see [6] for details) have performed a series of temperature-dependent 1 H measurements in order to characterize the hydrogen-bonding motifs involving the imidazopyridopyrimides bases. While they have been able to describe the formation of the amino-oxo–oxo-amino heterodimer and of the aminooxo and oxo-amino homodimers, rather complicated spectra have been obtained in the case of the diamino–dioxo heterodimer. Here, the HF/6-31G** geometries of the complexes 1 2, 3 4, 3 3, and 4 4 were utilized to predict the positions of the signals of the amino, amide, H2, and H9 protons. The calculated and experimental data are summarized in Table 3. First, the linear relationship was established between the computed absolute 1 H isotropic chemical shielding, r, of the above-mentioned nuclei in the 3 4, 3 3, and 4 4 dimers and the corresponding experimental chemical shifts, d, obtained at )60 °C. A very good correlation between the calculated and measured data was obtained, namely, r ¼ 1:193 d + 36.68 ppm, r2 ¼ 0.9761, rmsd ¼ 0.3985. Hence, this relationship was subsequently used to predict the corresponding 1 H spectrum of the 1 2 dimer. Based on these results, the assignment of the signals detected at )60 °C in CDCl3 solution (see supplementary information of [6]) is proposed (Table 3). Despite a number of approximations made during the theoretical treatment of the 1 H chemical shifts (see [12] for a discussion), the SOSDFPT-IGLO calculations provide fairly reliable estimates of the position of the peaks in experimental spectra. It is hoped that these results can be useful in interpreting potentially important measurements extended toward lower temperatures (which probably could be performed in, for example, dichloromethane). The intrinsic energies of the model imidazopyridopyrimide base pairs will now be put into perspective of thermodynamic measurements carried by Minakawa et al. [6] for the duplexes of oligodeoxynucleotides containing diamino–dioxo and amino-oxo–oxo-amino dimers. Minakawa et al. have shown experimentally that the hydrogen-bonding abilities of the imidazopyridopyrimide bases are essential for the thermal and thermodynamic stabilization of the duplexes of respective oligodeoxynucleotides. In particular, they have demonstrated that consecutive incorporation of diamino–dioxo and amino-oxo–oxo-amino pairs increases the melting temperatures of the duplexes by 4.2 and 18.2 °C, respectively, with respect to the GC, and by 10.1 and
512
J. Czernek / Chemical Physics Letters 392 (2004) 508–513
Table 3 The experimental and theoretical 1 H NMR data (ppm) of different complexesa Complex
Atom
r (calcd)
r (param)
d (exptl)
34
H5 (4) H6 (3) H4a (3) H4a (4) H9 (3) H9 (4) H2 H4b (3) H4b (4)
19.721 19.823 21.999 22.411 22.488 22.916 23.834 25.221 26.243
13.15 13.03 10.43 9.94 9.85 9.33 8.24 6.58 5.36
13.58 13.37 10.24 9.49 9.3 9.1 8.1 6.78 5.68
33
H6 H4a H9 H2 H4b
19.588 21.763 22.479 23.822 25.228
13.31 10.71 9.86 8.25 6.58
13.9 10.45 9.3 8.15 7.25
44
H5 H7a H9 H2 H7b
19.966 22.625 22.909 23.849 26.221
12.85 9.68 9.34 8.22 5.39
13.14 9.33 9.1 8.05 5.79
12
H5 H6 H4a H7a H9 (2) H9 (1) H2 (1) H2 (2) H4b H7b
18.251 18.394 22.059 22.592 22.621 22.778 23.827 23.845 25.263 26.312
14.90 14.73 10.36 9.72 9.69 9.50 8.25 8.23 6.53 5.28
14.8 14.8 10.4 9.4 9.4 9.1 8.1 8.1 6.4 5.5
a r (calcd) and r (param) denote the 1 H isotropic chemical shielding computed by the SOS-DPFT-IGLO method and obtained from the parametrization as described in the text, respectively, and d (exptl) is the measured value of the 1 H chemical shift.
23.5 °C, respectively, with respect to the AT Watson– Crick base pairs. Obviously, these results are contextdependendent and cannot be directly compared to the data obtained for the reduced systems in the gas phase, which were described in this Letter. Nevertheless, the intrinsic molecular interactions in imidazopyridopyrimide base pairs as revealed by the quantum chemical calculations agree qualitatively with the experimental findings and fully support the notion of the stable pairing motif consisting of four hydrogen bonds in the diamino–dioxo and amino-oxo–oxo-amino dimers.
4. Conclusions A high-level quantum chemical investigation of the intrinsic molecular interactions in imidazopyridopyrimidine base pairs was performed. A significant stabilization of the complexes containing four hydrogen bonds is predicted by all the methods employed. In particular, at the B3LYP/6-311G** level, the interaction energies of
146.7 and 102.5 kJ/mol have been found for the diamino–dioxo and amino-oxo–oxo-amino complexes, respectively, while the standard free energy changes at 298 K were established to be 63.3 and 27.9 kJ/mol accordingly. This qualitatively agrees with the results of thermodynamic measurements carried by Minakawa et al. [6] for the duplexes containing these dimers. The 1 H NMR chemical shifts of the complexes were predicted by the SOS-DFPT-IGLO method. Based on these calculations, the assignment of the low-temperature spectrum of the diamino–dioxo pair was proposed.
Acknowledgements This research has been supported by the Academy of Sciences of the Czech Republic, Grants AVOZ 4050913 and KJB4050311. Time allocation in the Czech Academic Supercomputer Centre and in the Mississippi Center for Supercomputing Research is gratefully acknowledged.
References [1] P. Belmont, J.-F. Constant, M. Demeunynck, Chem. Soc. Rev. 30 (2001) 70. [2] S. Scheiner, Hydrogen Bonding: A Theoretical Perspective, Oxford University Press, New York, 1997. [3] J.D. Watson, F.H.C. Crick, Nature 171 (1953) 737. [4] K. Hoogsteen, Acta Crystallogr. 16 (1963) 907. [5] L.M. Wilhelmsson, A. Holmen, P. Lincoln, P.E. Nielsen, B. Norden, J. Am. Chem. Soc. 123 (2001) 2434. [6] N. Minakawa, N. Kojima, S. Hikishima, T. Sasaki, A. Kiyosue, N. Atsumi, Y. Ueno, A. Matsuda, J. Am. Chem. Soc. 125 (2003) 9970. [7] J. Sponer, J. Leszczynski, P. Hobza, Biopolymers 61 (2002) 3. [8] D. Reha, M. Kabelac, F. Ryjacek, J. Sponer, J.E. Sponer, M. Elstner, S. Suhai, P. Hobza, J. Am. Chem. Soc. 124 (2002) 3366. ckova, J. Leszczynski, J. Sponer, [9] N. Gresh, J.E. Sponer, N. Spa J. Phys. Chem. B 107 (2003) 8669. [10] P. Hobza, J. Sponer, Chem. Rev. 99 (1999) 3247. [11] T. Helgaker, M. Jaszu nski, K. Ruud, Chem. Rev. 99 (1999) 293. [12] J. Czernek, V. Sklenar, J. Phys. Chem. A 103 (1999) 4089. [13] R. Fiala, J. Czernek, V. Sklenar, J. Biomol. NMR 16 (2000) 291. [14] D. Sitkoff, D.A. Case, Prog. NMR Spectrosc. 32 (1998) 165. [15] J. Czernek, R. Fiala, V. Sklenar, J. Magn. Res. 145 (2000) 142. [16] J. Czernek, J. Phys. Chem. A 105 (2001) 1357. [17] J. Sponer, J. Leszczynski, P. Hobza, J. Phys. Chem. 100 (1996) 1965. [18] O. Lukin, J. Leszczynski, J. Phys. Chem. A 106 (2002) 6775. [19] The optimized geometries, their total energies, the geometrical parameters of the hydrogen bonds, and the chemical shielding tensors of the structures investigated can be obtained from the author upon request. [20] A.D. Becke, J. Chem. Phys. 98 (1993) 5648. [21] C. Lee, W. Yang, R. Parr, Phys. Rev. B 37 (1998) 785. [22] S.F. Boys, F. Bernardi, Mol. Phys. 19 (1970) 553. [23] M.J. Frish et al., GA U S S I A N 98, Revision A7, Gaussian, Inc., Pittsburg, PA, 1998.
J. Czernek / Chemical Physics Letters 392 (2004) 508–513 [24] P. Hobza, J. Sponer, Chem. Phys. Lett. 261 (1996) 379. [25] V.G. Malkin, O.L. Malkina, M.E. Casida, D.R. Salahub, J. Am. Chem. Soc. 116 (1994) 5898. [26] D.R. Salahub, R. Fournier, P. Mlynarski, A. Papai, A. St-Amant, J. Uskio, in: J. Labanowski, J.W. Andzelm (Eds.), Density Functional Methods in Chemistry, Springer, New York, 1991, p. 77. [27] V.G. Malkin, O.L. Malkina, D.R. Salahub, MASTER-CS Program, Universite de Montreal, Montreal, 1994. [28] W. Kutzelnigg, U. Fleischer, M. Schindler (Eds.), NMR Basic Principles and Progress, vol. 23, Springer-Verlag, Berlin, 1990, p. 165. [29] J.P. Perdew, Y. Wang, Phys. Rev. B 45 (1992) 13244.
513
[30] J.P. Perdew, J.A. Chevary, S.H. Vosko, K.A. Jackson, M.R. Pederson, D.J. Singh, C. Fiolhais, Phys. Rev. B 46 (1992) 6671. [31] V.G. Malkin, O.L. Malkina, L.A. Eriksson, D.R. Salahub, in: J.M. Seminario, P. Politzer (Eds.), Theoretical and Computational Chemistry, vol. 2, Elsevier, Amsterdam, 1995, p. 273. [32] S.F. Boys, in: P.-O. L€ owdin (Ed.), Quantum Theory of Atoms, Molecules, and the Solid State, Academic Press, New York, 1969, p. 263. [33] J. Sponer, P. Hobza, Chem. Phys. Lett. 267 (1997) 263. [34] S.-I. Kawahara, T. Uchimaru, K. Tatra, M. Sekine, J. Phys. Chem. A 106 (2002) 3207. [35] B.G. Johnson, P.M.W. Gill, J.A. Pople, J. Chem. Phys. 98 (1995) 5612.
Imidazopyridopyrimidine base pairing motifs consisting of four hydrogen bonds: a quantum chemical study Jirı Czernek
*
Institute of Macromolecular Chemistry, Academy of Sciences of the Czech Republic, Heyrovsky Square 2, 162 06 Praha, The Czech Republic Received 28 January 2004; in final form 16 April 2004 Available online 19 June 2004
Abstract Following a recent experimental report on the synthesis and thermal stability of imidazopyridopyrimidine-containing structures with the ability to form four hydrogen bonds, the geometries and relative stabilities of the pairs of model bases have been studied using the HF, MP2 and B3LYP approaches. The 1 H NMR chemical shifts of the base pairs have been estimated by the SOS-DFPTIGLO method and employed in the analysis of the experimental proton spectra. Ó 2004 Elsevier B.V. All rights reserved.
1. Introduction DNA is the structural support of genetic material and therefore, the key factor in many vital cellular processes. While double-stranded right-handed helix is a regular conformation adopted by DNA in cells, there is an increasing number of results pointing to the biological importance of alternative structures such as bulges, hairpins, cruciforms, branched junctions, quadruplexes, and left-handed Z conformations [1]. Hydrogen bonding (H-bonding) [2] is of paramount importance in stabilizing helical (through the Watson–Crick base pairing [3]) and also the other DNA structures, which involve alternative pairing schemes such as Hoogsteen type interactions [4]. Consequently, a considerable attention has been paid to the design, synthesis and structural analysis of short oligodeoxynucleotides (ODNs) featuring different H-bonding motifs [5]. Most recently, two sets of new base pair combinations consisting of four hydrogen bonds have been proposed (see Fig. 1) and incorporated into ODNs [6]. Importantly, the tricyclic skeletons adopted by Minakawa et al., [6] i.e., imidazo[50 ; 40 :4,5]pyrido[2,3-d] pyrimidine nucleosides, possess the ability to form both Watson–Crick and Hoogsteen types of H-bonding and thus are potentially *
Tel.: +420-296809290; Fax: +420-296809410. E-mail address: [email protected] (J. Czernek).
0009-2614/$ - see front matter Ó 2004 Elsevier B.V. All rights reserved. doi:10.1016/j.cplett.2004.06.007
useful in stabilizing both the double-stranded helical and alternative DNA structures. These novel pairing motifs have been characterized using thermal stabilities of several duplexes of ODNs containing different number of the tricyclic nucleosides. It has been shown that in some cases the thermal stability of imidazopyridopyrimide pairs surpasses that of their Watson–Crick counterparts. Presently, high-level ab initio computational techniques are able to reliably account for the intrinsic interactions in nucleic acids, i.e., to quantitatively describe the nature and energetics of various interactions in (sizeable fragments of) nucleic acids without the presence of external effects such as solvent, crystal packing and others (see [7] for a recent review). The knowledge of the intrinsic molecular interactions is indispensable in understanding the role played by respective structural elements, including intercalators [8] and metal cations [9], in DNA structure and function. Thanks to a significant benchmarking effort, useful guidelines for correct and computationally efficient computational studies of nucleic acids have been devised [10]. As these approaches also provide for accurate geometrical characteristics of the studied structures, an advantage can be taken from recent advances in ab initio methods for calculating NMR (nucleic magnetic resonance) parameters [11] to predict and interpret the spectra of nucleosides [12,13] and base pairs [14–16].
J. Czernek / Chemical Physics Letters 392 (2004) 508–513
(a)
H4b
H2
N
N
O
H4a
H2
N
HN
N
H5
N
NH
H9
N
H6
N
H9 N
N N
O
H7a
H7b
diamino (1)
dioxo (2)
H4b
(b)
H2
N
N HN
N
H9
N
O
H4a
H5
H6
H2
N
N
NH
N
H9 N
N O
H7a
N H7b
amino-oxo (3)
oxo-amino (4)
Fig. 1. Base pairing in diamino–dioxo (1 2,a) and amino-oxo–oxoamino (3 4,b) complexes. The numbering of protons for which the chemical shielding has been evaluated is shown.
The main objective of this Letter, is the quantum chemical description of the structure and stability of the diamino–dioxo and amino-oxo–oxo-amino heterodimers of imidazopyridopyrimide base pairs (Fig. 1). The results are compared with the computational data obtained at the same, sufficiently high levels of theory for guanine–cytosine (GC) and adenine–thymine (AT) Watson–Crick base pairs [17], and for selected dimers containing four hydrogen bonds [18], and are discussed in the context of the thermal stabilities found by Minakawa et al. [6]. Moreover, to aid in the interpretation of the 1 H NMR spectra of different pairing motifs, the amino-oxo and oxo-amino homodimers have been studied in addition to the above-mentioned heterodimers, and the chemical shifts of all the base pairs have been computed using a density functional theory (DFT)based method.
2. Calculations The structures of four dimers formed by imidazo[50 ; 40 :4,5]pyrido[2,3-d]pyrimidine 1 (denoted as ‘diamino’ in [6]), imidazo[50 ; 40 :4,5]pyrido[2,3-d]pyrimi-
509
dine-4,7(5H,6H)-dione 2 (‘dioxo’), imidazo[50 ; 40 :4,5]pyrido[2,3-d]pyrimidin-7(6H)-one 3 (‘amino-oxo’), and imidazo[50 ; 40 :4,5]pyrido[2,3-d]pyrimidin-4(5H)-one 4 (‘oxo-amino’) were studied [19], namely, the heterodimers 1 2 and 3 4, and the homodimers 3 3 and 4 4 (see Fig. 1). Thus, the molecular geometries of these dimers and of the monomers of 1–4 were optimized using the Hartree–Fock (HF) theory and the 6-31G** basis set of atomic orbitals, and the density functional theory with B3LYP (Becke’s three parameter exchange [20] and Lee, Yang, Parr correlation [21]) functionals and the 6-311G** basis set. The Cs symmetry was imposed on the monomers and heterodimers, while the C2h symmetry was exploited in the calculations of the homodimers. To verify the character of the optimized structures and to obtain the data for thermodynamic computations (see below), calculations of the harmonic vibrational frequencies were subsequently performed at the same level as the geometrical optimizations. The interaction energies, DE, were calculated using the variational supermolecular approach as the energy difference between the complex and the isolated components. Namely, for the structures optimized at the HF/6-31G** level, the second-order Møller–Plesset (MP2) energies were obtained in the frozen core approximation using the standard 6-31G basis set augmented with one set of diffuse d-polarization functions with the exponent of 0.25, which were placed on the second-row atoms. This methodology will be referred to as MP2/6-31G*(0.25)// HF/6-31G**. To assess the importance of correlation energy for the stabilization of imidazopyridopyrimide pairs, the HF/6-31G*(0.25) interaction energies were evaluated as well. Moreover, for the structures optimized using the B3LYP/6-311G** approach, the DE’s were computed at the same level of theory. All interaction energies were corrected for the basis set superposition error by means of the counterpoise method of Boys and Bernardi [22]. In addition, for the B3LYP/ 6-311G** data the deformation energy, EDEF , was computed by subtracting the energy of the optimized isolated monomers from the energies of the monomers in the dimer geometry [10]. All the above-mentioned computations were performed using the default algorithms and settings of GA U S S I A N 98 [23]. Employing the interaction energies and the unscaled values of harmonic vibrational frequencies, the B3LYP/6-311G** enthalpies and free energies of the complexes at 298 K were estimated in the rigid rotor-harmonic oscillator-ideal gas approximation [24]. The 1 H NMR chemical shielding tensors were obtained using the SOS-DFPT-IGLO [25] method as implemented in the deMon-MASTER-CS code [26,27]. The SOS-DFPT-IGLO strategy combines the sum-overstates density functional (Rayleigh–Schr€ odinger) perturbation theory with the individual gauges for localized (molecular) orbitals [28]. The IGLO-III basis set of
510
J. Czernek / Chemical Physics Letters 392 (2004) 508–513
Kutzelnigg et al. [28] (roughly of quadruple-n quality with two sets of polarization functions), the Perdew– Wang [29,30] exchange-correlation potential, the FINE angular integration grid with 64 radial shells [31], and the approximation Loc. 1 SOS-DFPT [25] were applied. The molecular orbitals were localized by the method of Boys [32].
3. Results and discussion Table 1 summarizes the counterpoise–corrected interaction energies as provided by the well-tested [10] MP2/6-31G*(0.25)//HF/6-31G** strategy for the 1 2, 3 4, 3 3, and 4 4 dimers. DE’s of an important type of structures which have been studied using the abovementioned approach, i.e., the GC and AT Watson– Crick base pairs [17], are shown for comparison purposes. The data obtained for GC and AT have been shown, due to a compensation of errors, to lie close to the results of the CCSD(T)/aug-cc-pVDZ calculations [33]. Thus, the MP2/6-31G*(0.25)//HF/6-31G** approach is expected to yield reliable characteristics also of the imidazopyridopyrimide bases. Data in Table 1 clearly indicate a significant stabilizing interaction in all the dimers investigated, with DE’s ranging from 91.3 (in the case of 3 3) to 136.0 kJ/mol (1 2). These values are comparable to the complexation energy of the GC pair, but they are ca. two times higher than DE of the AT Watson–Crick base pair (see Table 1). It has been demonstrated that the stabilization of hydrogen-bonded bases of nucleic acids [10], and of similar types of hydrogen-bonding systems [34], is predominantly electrostatic in origin. Hence, quite reliable characteristics of this type of structures can be obtained already at the HF level with medium-sized polarized basis sets [33]. This holds also for four imidazopyridopyrimide base pairs studied here: the order of the HF/6-31G*(0.25) interaction energies is identical to that found by the MP2 method. The correlation contribution to the stabiliza-
Table 1 The negative of total MP2 interaction energy, DEMP2 , and its HF, DEHF , and electron correlation, DECOR , components (kJ/mol), calculateda for different complexes Complex
DEMP2
DEHF
DECOR
12 34 33 44 GC AT
136.0 98.0 105.5 91.3 108.4 52.1
117.8 78.7 85.0 72.9 103.3 40.7
18.1 19.3 20.5 18.4 5.0 11.3
See the text for details. a The results for the GC and AT base pairs are taken from Table 1 of [17].
tion energy, i.e., the difference between the respective MP2 and HF data, constitutes 20% of the DEMP2 in the case of 3 4, 3 3 and 4 4 pairs, and 13% in the case of the most stable dimer 1 2. This contribution amounts to 5% and 22% for the GC and AT pairs, respectively (Table 1). It is noted that all vibrational frequencies obtained for the HF/6-31G** geometries have been found to be real except the one reflecting inversion of the amino group [10] in the position 7 of the diamino monomer 1. Let us now turn to the discussion of the B3LYP/ 6-311G** results. This approach, which now can be routinely applied to quite sizeable systems, is known to be highly accurate. For example, it has been shown to be capable of providing smaller absolute deviations for the bond lengths and angles than the QCISD/6-31G* method [35], and successfully applied in a recent comparative study of the geometries and relative stabilities of 17 multiply hydrogen-bonded complexes [18]. Consequently, the B3LYP/6-311G** technique was adopted here to study the geometry, energetics, and thermodynamics of the imidazopyridopyrimide base pairs. As was the case of the HF/6-31G** geometries, all B3LYP/ 6-311G** vibrational frequencies are real except of the one describing inversion of the amino group in the position 7 of the compound 1. The data relevant to the thermodynamics of imidazopyridopyrimide dimers formation are summarized in Table 2 together with the results for the most and the least stable complexes containing four hydrogen bonds, taken from the set recently studied by Lukin and Leszczynski, [18] i.e., the structures denoted as 11 and 14, respectively (see also Fig. 2). In accord with the MP2/6-31G*(0.25)//HF/6-31G**computed interaction energies, a significant stabilization of imidazopyridopyrimide dimers is predicted on the basis of the B3LYP/6-311G** calculations. These two approaches also agree in the ordering of the DE’s for respective pairs, with the 3 3 being the least and the 1 2 the most stable dimers (their DE’s amount to 96.1 and Table 2 The negative of the DFT-computed: interaction energy, DEDFT , deformation energy, EDEF , and standard molar enthalpy, DH 0 298 , entropy, DS 0 298 , and Gibbs energy, DG0 298 , changes at 298 K calculateda for different complexes (kJ/mol)b Complex
DEDFT
EDEF
DH 0 298
DS 0 298
DG0 298
12 34 33 44 11 14
146.7 102.5 110.0 96.1 196.4 58.1
)17.5 )9.9 )11.2 )8.9 NA NA
119.5 82.9 91.6 75.3 191.3 69.7
188.2 184.4 201.6 173.5 249.3 213.8
63.3 27.9 31.5 23.6 117.0 6.0
See the text for details. The results for the complexes denoted as 11 and 14 are taken from Table 2 of [18]. b DS 0 298 is given in J/(mol K). a
J. Czernek / Chemical Physics Letters 392 (2004) 508–513
H O
N O
N
H
H
N
N
N
N
H
H
N
O
CH3
H3C O
N H
11 O H
H N H O
N H N CH3
N
N
H
H
N
O NH
14 Fig. 2. The structure of the complexes denoted as 11 and 14.
146.7 kJ/mol accordingly). In agreement with previous findings [10] is the fact that the larger the stabilization, the higher is the (always repulsive) deformation energy of the bases (see Table 2). An inclusion of the zero-point vibration energy and temperature-dependent enthalpy contributions moderately reduces the stabilization (as expressed by the values of DH 00 and DH 0298 accordingly). The entropic term is fairly uniform (16% variation within the series), but significant, constituting up to 69% of the standard molar enthalpy change at 298 K. However, its inclusion does not alter the stability order of imidazopyridopyrimide base pairs as expressed by their DG0298 ’s (Table 2). Clearly, apart from the solvent effects, there are many structural factors affecting the stability of hydrogenbonded systems, including the type of donor–acceptor sequences, specific substituents, and conformational properties. Nevertheless, it can be of some interest to compare the results obtained for imidazopyridopyrimide base pairs with those calculated at the B3LYP/ 6-311G** level for some complexes containing four hydrogen bonds [18]. The energetic and thermodynamic data describing the most and the least stable assemblies studied by Lukin and Leszczynski [18] generally encompass the values found for imidazopyridopyrimide dimers. However, the exception is the (destabilizing) entropic term, which is for imidazopyridopyrimides bases smaller than for both the dimers 11 and 14 (see
511
Fig. 2 and Table 2). In the case of the least stable structure 14 from [18], the enthalpy stabilization is almost cancelled out by the entropy contribution, leading to the standard Gibbs energy change at 298 K of only 6.0 kJ/mol. On the other hand, the structure 11 is significantly more stable than imidazopyridopyrimide base pairs, partly due to the presence of two intramolecular hydrogens bonds [18]. Minakawa et al. (see [6] for details) have performed a series of temperature-dependent 1 H measurements in order to characterize the hydrogen-bonding motifs involving the imidazopyridopyrimides bases. While they have been able to describe the formation of the amino-oxo–oxo-amino heterodimer and of the aminooxo and oxo-amino homodimers, rather complicated spectra have been obtained in the case of the diamino–dioxo heterodimer. Here, the HF/6-31G** geometries of the complexes 1 2, 3 4, 3 3, and 4 4 were utilized to predict the positions of the signals of the amino, amide, H2, and H9 protons. The calculated and experimental data are summarized in Table 3. First, the linear relationship was established between the computed absolute 1 H isotropic chemical shielding, r, of the above-mentioned nuclei in the 3 4, 3 3, and 4 4 dimers and the corresponding experimental chemical shifts, d, obtained at )60 °C. A very good correlation between the calculated and measured data was obtained, namely, r ¼ 1:193 d + 36.68 ppm, r2 ¼ 0.9761, rmsd ¼ 0.3985. Hence, this relationship was subsequently used to predict the corresponding 1 H spectrum of the 1 2 dimer. Based on these results, the assignment of the signals detected at )60 °C in CDCl3 solution (see supplementary information of [6]) is proposed (Table 3). Despite a number of approximations made during the theoretical treatment of the 1 H chemical shifts (see [12] for a discussion), the SOSDFPT-IGLO calculations provide fairly reliable estimates of the position of the peaks in experimental spectra. It is hoped that these results can be useful in interpreting potentially important measurements extended toward lower temperatures (which probably could be performed in, for example, dichloromethane). The intrinsic energies of the model imidazopyridopyrimide base pairs will now be put into perspective of thermodynamic measurements carried by Minakawa et al. [6] for the duplexes of oligodeoxynucleotides containing diamino–dioxo and amino-oxo–oxo-amino dimers. Minakawa et al. have shown experimentally that the hydrogen-bonding abilities of the imidazopyridopyrimide bases are essential for the thermal and thermodynamic stabilization of the duplexes of respective oligodeoxynucleotides. In particular, they have demonstrated that consecutive incorporation of diamino–dioxo and amino-oxo–oxo-amino pairs increases the melting temperatures of the duplexes by 4.2 and 18.2 °C, respectively, with respect to the GC, and by 10.1 and
512
J. Czernek / Chemical Physics Letters 392 (2004) 508–513
Table 3 The experimental and theoretical 1 H NMR data (ppm) of different complexesa Complex
Atom
r (calcd)
r (param)
d (exptl)
34
H5 (4) H6 (3) H4a (3) H4a (4) H9 (3) H9 (4) H2 H4b (3) H4b (4)
19.721 19.823 21.999 22.411 22.488 22.916 23.834 25.221 26.243
13.15 13.03 10.43 9.94 9.85 9.33 8.24 6.58 5.36
13.58 13.37 10.24 9.49 9.3 9.1 8.1 6.78 5.68
33
H6 H4a H9 H2 H4b
19.588 21.763 22.479 23.822 25.228
13.31 10.71 9.86 8.25 6.58
13.9 10.45 9.3 8.15 7.25
44
H5 H7a H9 H2 H7b
19.966 22.625 22.909 23.849 26.221
12.85 9.68 9.34 8.22 5.39
13.14 9.33 9.1 8.05 5.79
12
H5 H6 H4a H7a H9 (2) H9 (1) H2 (1) H2 (2) H4b H7b
18.251 18.394 22.059 22.592 22.621 22.778 23.827 23.845 25.263 26.312
14.90 14.73 10.36 9.72 9.69 9.50 8.25 8.23 6.53 5.28
14.8 14.8 10.4 9.4 9.4 9.1 8.1 8.1 6.4 5.5
a r (calcd) and r (param) denote the 1 H isotropic chemical shielding computed by the SOS-DPFT-IGLO method and obtained from the parametrization as described in the text, respectively, and d (exptl) is the measured value of the 1 H chemical shift.
23.5 °C, respectively, with respect to the AT Watson– Crick base pairs. Obviously, these results are contextdependendent and cannot be directly compared to the data obtained for the reduced systems in the gas phase, which were described in this Letter. Nevertheless, the intrinsic molecular interactions in imidazopyridopyrimide base pairs as revealed by the quantum chemical calculations agree qualitatively with the experimental findings and fully support the notion of the stable pairing motif consisting of four hydrogen bonds in the diamino–dioxo and amino-oxo–oxo-amino dimers.
4. Conclusions A high-level quantum chemical investigation of the intrinsic molecular interactions in imidazopyridopyrimidine base pairs was performed. A significant stabilization of the complexes containing four hydrogen bonds is predicted by all the methods employed. In particular, at the B3LYP/6-311G** level, the interaction energies of
146.7 and 102.5 kJ/mol have been found for the diamino–dioxo and amino-oxo–oxo-amino complexes, respectively, while the standard free energy changes at 298 K were established to be 63.3 and 27.9 kJ/mol accordingly. This qualitatively agrees with the results of thermodynamic measurements carried by Minakawa et al. [6] for the duplexes containing these dimers. The 1 H NMR chemical shifts of the complexes were predicted by the SOS-DFPT-IGLO method. Based on these calculations, the assignment of the low-temperature spectrum of the diamino–dioxo pair was proposed.
Acknowledgements This research has been supported by the Academy of Sciences of the Czech Republic, Grants AVOZ 4050913 and KJB4050311. Time allocation in the Czech Academic Supercomputer Centre and in the Mississippi Center for Supercomputing Research is gratefully acknowledged.
References [1] P. Belmont, J.-F. Constant, M. Demeunynck, Chem. Soc. Rev. 30 (2001) 70. [2] S. Scheiner, Hydrogen Bonding: A Theoretical Perspective, Oxford University Press, New York, 1997. [3] J.D. Watson, F.H.C. Crick, Nature 171 (1953) 737. [4] K. Hoogsteen, Acta Crystallogr. 16 (1963) 907. [5] L.M. Wilhelmsson, A. Holmen, P. Lincoln, P.E. Nielsen, B. Norden, J. Am. Chem. Soc. 123 (2001) 2434. [6] N. Minakawa, N. Kojima, S. Hikishima, T. Sasaki, A. Kiyosue, N. Atsumi, Y. Ueno, A. Matsuda, J. Am. Chem. Soc. 125 (2003) 9970. [7] J. Sponer, J. Leszczynski, P. Hobza, Biopolymers 61 (2002) 3. [8] D. Reha, M. Kabelac, F. Ryjacek, J. Sponer, J.E. Sponer, M. Elstner, S. Suhai, P. Hobza, J. Am. Chem. Soc. 124 (2002) 3366. ckova, J. Leszczynski, J. Sponer, [9] N. Gresh, J.E. Sponer, N. Spa J. Phys. Chem. B 107 (2003) 8669. [10] P. Hobza, J. Sponer, Chem. Rev. 99 (1999) 3247. [11] T. Helgaker, M. Jaszu nski, K. Ruud, Chem. Rev. 99 (1999) 293. [12] J. Czernek, V. Sklenar, J. Phys. Chem. A 103 (1999) 4089. [13] R. Fiala, J. Czernek, V. Sklenar, J. Biomol. NMR 16 (2000) 291. [14] D. Sitkoff, D.A. Case, Prog. NMR Spectrosc. 32 (1998) 165. [15] J. Czernek, R. Fiala, V. Sklenar, J. Magn. Res. 145 (2000) 142. [16] J. Czernek, J. Phys. Chem. A 105 (2001) 1357. [17] J. Sponer, J. Leszczynski, P. Hobza, J. Phys. Chem. 100 (1996) 1965. [18] O. Lukin, J. Leszczynski, J. Phys. Chem. A 106 (2002) 6775. [19] The optimized geometries, their total energies, the geometrical parameters of the hydrogen bonds, and the chemical shielding tensors of the structures investigated can be obtained from the author upon request. [20] A.D. Becke, J. Chem. Phys. 98 (1993) 5648. [21] C. Lee, W. Yang, R. Parr, Phys. Rev. B 37 (1998) 785. [22] S.F. Boys, F. Bernardi, Mol. Phys. 19 (1970) 553. [23] M.J. Frish et al., GA U S S I A N 98, Revision A7, Gaussian, Inc., Pittsburg, PA, 1998.
J. Czernek / Chemical Physics Letters 392 (2004) 508–513 [24] P. Hobza, J. Sponer, Chem. Phys. Lett. 261 (1996) 379. [25] V.G. Malkin, O.L. Malkina, M.E. Casida, D.R. Salahub, J. Am. Chem. Soc. 116 (1994) 5898. [26] D.R. Salahub, R. Fournier, P. Mlynarski, A. Papai, A. St-Amant, J. Uskio, in: J. Labanowski, J.W. Andzelm (Eds.), Density Functional Methods in Chemistry, Springer, New York, 1991, p. 77. [27] V.G. Malkin, O.L. Malkina, D.R. Salahub, MASTER-CS Program, Universite de Montreal, Montreal, 1994. [28] W. Kutzelnigg, U. Fleischer, M. Schindler (Eds.), NMR Basic Principles and Progress, vol. 23, Springer-Verlag, Berlin, 1990, p. 165. [29] J.P. Perdew, Y. Wang, Phys. Rev. B 45 (1992) 13244.
513
[30] J.P. Perdew, J.A. Chevary, S.H. Vosko, K.A. Jackson, M.R. Pederson, D.J. Singh, C. Fiolhais, Phys. Rev. B 46 (1992) 6671. [31] V.G. Malkin, O.L. Malkina, L.A. Eriksson, D.R. Salahub, in: J.M. Seminario, P. Politzer (Eds.), Theoretical and Computational Chemistry, vol. 2, Elsevier, Amsterdam, 1995, p. 273. [32] S.F. Boys, in: P.-O. L€ owdin (Ed.), Quantum Theory of Atoms, Molecules, and the Solid State, Academic Press, New York, 1969, p. 263. [33] J. Sponer, P. Hobza, Chem. Phys. Lett. 267 (1997) 263. [34] S.-I. Kawahara, T. Uchimaru, K. Tatra, M. Sekine, J. Phys. Chem. A 106 (2002) 3207. [35] B.G. Johnson, P.M.W. Gill, J.A. Pople, J. Chem. Phys. 98 (1995) 5612.