J. Phys. Chem. B 2009, 113, 5653–5656
5653
Solvent Reorganization Energy of Hole Transfer in DNA Toma´sˇ Kubarˇ† and Marcus Elstner*,†,‡ Institute of Physical Chemistry, Technische UniVersita¨t Braunschweig, D-38106 Braunschweig, Germany, and Department of Molecular Biophysics, German Cancer Research Center, D-69115 Heidelberg, Germany ReceiVed: December 21, 2008
We present an application of a molecular-dynamics-based scheme to evaluate the solvent reorganization energy of hole transfer in DNA. The obtained parameters can be used for simulations of hole transfer in DNA by means of Marcus’ theory. Also, we perform an analysis of the reorganization energies, including the case of transfer of a delocalized hole. I. Introduction In the past decade, the long-distance charge transfer (CT) in DNA stands in the center of interest of both experimental and computational research efforts in chemistry and physics.1 In computational studies (as in refs 2-5), the kinetics of CT in DNA is often described by means of Marcus’ theory,6,7 considering the nucleobases (and sometimes DNA ligands) as charge donors and acceptors and using the well-known equations for the rate of transfer
[
(∆G0 + λ)2 k ) κν exp 4λkBT
]
(1)
in the adiabatic regime (κ ) electronic transmission coefficient, ν ) effective frequency for nuclear motion along the reaction coordinate, kB ) Boltzmann constant, T ) temperature) or
k)
[
(∆G0 + λ)2 1 2π |HDA | 2 exp p 4λkBT √4πλkBT
]
(2)
under the assumption of weak electronic coupling (p ) Dirac constant, HDA ) electronic coupling). The activation-energy term contains the free energy change of charge transfer ∆G0 and the reorganization energy (RE) λ. The latter represents the amount of energy required to reorganize the donor and acceptor molecules as well as their environment (solvent and backbone in the case of DNA) from their initial configuration to that after CT has taken place. If the donor and acceptor of excess charge are chemically identical (e. g., both are guanines), then ∆G0 vanishes. In such a case, λ becomes the crucial quantity codetermining the CT rate. It may be decomposed into two contributions λ ) λi + λs: the internal RE λi, brought about by the change of molecular structure of charge donor and acceptor, and the solvent RE λs, which corresponds to the reorganization of the surrounding medium. λi is independent of the donoracceptor distance and may be easily evaluated using quantumchemical methods.8 On the other hand, the distance dependence of λs determines (together with HDA) how the rate of CT * Electronic address:
[email protected]. Phone: +49-531-3915347. Fax: +49-531-3915832. † Technische Universita¨t Braunschweig. ‡ German Cancer Research Center.
Figure 1. Reorganization energy λ.
decreases with distance, and this is more difficult to estimate, both experimentally and computationally. Several computational approaches were developed in order to estimate λs of CT in biomolecular systems such as DNA. First, the original Marcus picture of a continuum solvent was augmented by partitioning of the environment to regions of different permittivities, taking into account different electric properties of DNA bases, backbone, and surrounding water.9-12 Here, a correct choice of these permittivity values is rather difficult; the different sets of permittivities yield similar distance dependences of λs, yet their absolute magnitudes remain quite different. Second, explicit molecular-mechanics (MM) models were used to describe the interaction of solvent with the charge donor and acceptor represented by simple spherical particles, as well as with peptide complexes. λs was then evaluated for the charge separation and recombination reactions.13-16 These approaches avoid the issue of uncertainty of solvent and bridge permittivities. However, to our knowledge, there is no application of such a framework to CT in DNA. Our aim in this article is to provide a reliable estimate of the solvent reorganization energy of hole transfer in DNA, considering the explicit molecular environment of DNA by means of an MM force field. II. Methodology λsof the transfer from state R to β (differing by the placement of the hole) will be evaluated directly from the definition, as the increase of potential energy calculated with the Hamiltonian of state R when reorganizing the structure of the environment to that corresponding to state β: λs ) 〈VR(r)〉β - 〈VR(r)〉R (the subscript to V indicates the Hamiltonian used to calculate the energy, whereas the subscript to parentheses 〈〉 identifies the Hamiltonian used to create the statistical ensemble; cf. also Figure 1).
10.1021/jp901888r CCC: $40.75 2009 American Chemical Society Published on Web 03/30/2009
5654 J. Phys. Chem. B, Vol. 113, No. 16, 2009
Kubarˇ and Elstner
In this framework, we take advantage of molecular-dynamics (MD) simulation to provide sufficient sampling in the evaluation of energy. So, we obtain the energy 〈VR(r)〉R as the mean potential energy in the simulation of state R, whereas the energy 〈VR(r)〉β is calculated in a more complex procedure: an MD simulation is run with the Hamiltonian Vβ, and the desired energy is evaluated with the Hamiltonian VR. In practice, we first perform a normal MD simulation of the initial state of hole transfer and record the potential energy. Then, we perform a simulation of the final state and use this trajectory to re-evaluate the potential energy with the Hamiltonian corresponding to the initial state. The simulation lengths accessible on mainstream desktop computers, on the order of tens nanoseconds, are fully sufficient to guarantee adequate sampling and thus negligible statistical errors. The Amber 9 suite17 with the parm99+parmBSC0 force field18,19 was used to generate the topology and MM parameters of DNA species. The MM atomic charges on guanine and adenine carrying the hole (or its part) were estimated as follows: First, the atomic charges of a neutral nucleobase {q0i } as well as its cation {q+ i } were determined using the RESP methodology20 on the HF/6-31G* level of theory, consistently with the parm99 parametrization. So, the map of a hole on the nucleobase was obtained:
∀atomsi: ∆qi ) qi+ - qi0 ;
∑ ∆qi ) 1
(3)
i
Finally, the charges {qMM i } on MM residues were modified by addition of the corresponding fraction of the hole charge Qhole, which equals 1 for a hole localized on the respective nucleobase and is smaller than 1 if the hole is spread over several bases:
∀atomsi: qiMM′ ) qiMM + Qhole∆qi
(4)
However approximative, this approach made us able to describe the nonuniform distribution of a hole on the nucleobases as well as the environment of nucleobases in atomic detail. The simulations as well as potential energy re-evaluations were performed with the Gromacs 4.0 package.21,22 Standard equilibration and simulation protocols were applied and the length of production runs was 50 ns. The mean values of potential energy yielded by these simulations were checked for convergence both by comparison with the results of shorter simulations and by testing the normality of resulting distributions. No convergence problems were detected, and the sampling is thus considered sufficient. III. Results and Discussion First, we studied the reorganization energy of hole transfer in double-stranded DNA of sequences G6, A6, and GTnG (n ) 0,..., 4, used in the experiments by Giese et al.),23 capped by two guanines on each end. The transfer of an integral hole was considered from the 5′-end nucleobase to every other (G6 and A6) or to the 3′-end guanine (GTnG). This way, we obtained the distance dependence of λs, cf. Figure 2. λs increases with the distance between the hole donor and the acceptor, in accordance with the dependence reported in previous studies, both experimental24,25 and computational.9-16,30 In a recent study, Steinbrecher et al. estimated the λs of CT between neighboring nucleobases.26 The reported values of 1.08
Figure 2. Distance dependence of reorganization energy for hole transfer in various DNA sequences.
and 1.15 eV for a pair of adenines and guanines, respectively, are in reasonable agreement with our results of 1.21 and 1.41 eV, respectively. The modest discrepancy may have been brought about by the different hole-charge distribution on the nucleobases, which was caused probably by the more approximative quantum-chemical treatment by Steinbrecher et al. (only the π orbitals of the nucleobases were considered). Also, we observe small deviations between the examined DNA sequences, which may be attributed to their different helical structures (B-like in the case of poly(A), A-like in the case of poly(G)) as well as to the probably different hydration patterns and strengths of guanine and adenine. In further applications, the described calculation shall be repeated for the DNA sequence of interest and the obtained values of λs may be readily cast in the Marcus equation (eq 1 or 2). Some proposed mechanisms of CT in DNA1 involve a hole extended over several nucleobases rather than confined to a single one. Therefore, it is of interest to examine λs for the transfer of nonintegral excess charge (a part of a hole) as well as for the transfer of such a delocalized hole. As for the former case, the hole was considered extended over two guanines in the center of a G6 hexamer. Then, λs was evaluated for the transfer of partial charge from one of the guanines to the other. The results (cf. Figure 3, top) show a roughly quadratic dependence of reorganization energy on the amount of charge transferred. Thus, a transfer of a large amount of charge requires a higher reorganization barrier to be overcome than if sequential transfer of smaller charges were considered. Also, the dependence is not completely symmetric, which means that the energetic penalty stemming from reorganization in the directions 5′f3′ and 5′r3′ in the same DNA sequence may differ. In our example, λs amounts to 1.31 eV in the former direction and 1.14 eV in the latter, for the transfer of an integral hole. The potential energy itself is of interest, too, as it may tell us which charge states will be more favorable due to the electrostatic interactions with the solvent as well as within DNA. The mean values of total MM energy in the simulation of a differently (de)localized hole are presented in Figure 3, bottom. The difference in energy can be interpreted in terms of different solvation energies of the considered states of the hole, which were reported earlier.27 The energy difference amounts to ca. 0.7 eV for a hole shared in a 1:1 ratio by two nucleobases, relative to a fully localized hole. Thus, the delocalization of the hole seems to be penalized by a less favorable solvation of the states where the hole is spread over two nucleobases, as compared to a hole fully localized on one base or the other, in
Solvent Reorganization Energy of Hole Transfer in DNA
J. Phys. Chem. B, Vol. 113, No. 16, 2009 5655 ence smaller RE and thus a smaller activation barrier. The energy balance of the hole delocalization is quite complex; one should consider at least the following contributions: The electronic effects favor delocalizationsKurnikov et al.27 reported a hole spread over two nucleobases to be 0.47 eV more stable than the localized one. However, the localized hole is more favorably hydrated, and thus the delocalized hole lies 0.7 eV higher in solvation energy than the localized one. Also, λi for the transfer of partial charge from a nucleobase to its neighbor applies; this was evaluated at around 0.1 eV.31 Then, assuming the additivity of these effects in a naive fashion, the activation energy for hole delocalization might be estimated at around 0.3 eV. Once this barrier is overcome, CT of the delocalized hole may commence. Careful QM/MM-based treatment would be required to examine this phenomenon in a more accurate way. IV. Conclusion
Figure 3. Top: Reorganization energy of transfer of partial charge between neighboring guanines. Curves labeled by initial state; 0.0 and 1.0 denote the entire hole localized on the 5′-base and on the 3′-base, respectively. Bottom: Total MM energy in the simulations with a hole, in dependence on the degree of delocalization (arbitrarily set to zero in state 0.0).
In summary, we have applied a computational scheme to evaluate the solvent reorganization energy of hole transfer in DNA. The calculated values come within the interval predicted previously, as does their distance dependence. Also, we studied the processes of formation and transfer of a delocalized hole. Although the delocalization of a hole is energetically unfavorable, the transfer of such a hole spread over a few nucleobases experiences a lower energetic barrier than if the hole is confined to a single nucleobase. Acknowledgment. This work was supported by the Deutsche Forschungsgemeinschaft, Project DFG-EL 206/5-1. References and Notes
Figure 4. Reorganization energy of transfer of a delocalized hole in poly(G) DNA, compared to that of a localized one.
accord with the conclusions of Voityuk.28 Obviously, the process of relocalization of a hole would have to pass over a barrier of λs first, which may be estimated at 0.4 eV on the basis of our calculation (Figure 3, top). The other casesthat of transfer of a delocalized holeswas studied in the poly(G) sequence as well. The hole was placed on two neighboring guanines so that each of them bore a half of the hole charge. Then, λs was evaluated for the transfer of the hole to other pairs of neighboring bases, in the same way as in previous cases. The calculated values of λs are presented in Figure 4 together with the previous results for the localized hole. With the delocalized hole, we obtain a roughly 0.9 eV smaller λs than with the one confined to a single guanine, for every donor-acceptor distance. Therefore, if we assume the existence and transfer of a hole spread over several nucleobases, then this transfer will experi-
(1) Schuster, G. B., Ed. Long-Range Charge Transfer in DNA; Topics in Current Chemistry Vol. 236-237; Springer: 2004. (2) Jakobsson, M.; Stafstro¨m, S. J. Chem. Phys. 2008, 129, 125102. (3) Rˇeha, D.; Barford, W.; Harris, S. Phys. Chem. Chem. Phys. 2008, 10, 5436. (4) Cramer, T.; Krapf, S.; Koslowski, T. J. Phys. Chem. B 2004, 108, 11812. (5) Senthilkumar, K.; Grozema, F. C.; Guerra, C. F.; Bickelhaupt, F. M.; Lewis, F. D.; Berlin, Y. A.; Ratner, M. A.; Siebbeles, L. D. A. J. Am. Chem. Soc. 2005, 127, 14894. (6) Marcus, R. A. J. Chem. Phys. 1956, 24, 966. (7) Marcus, R. A.; Sutin, N. Biochim. Biophys. Acta 1985, 811, 265. (8) Berashevich, J. A.; Chakraborty, T. Chem. Phys. Lett. 2007, 446, 159. (9) Tavernier, H. L.; Fayer, M. D. J. Phys. Chem. B 2000, 104, 11541. (10) Tong, G. S. M.; Kurnikov, I. V.; Beratan, D. N. J. Phys. Chem. B 2002, 106, 2381. (11) Siriwong, K.; Voityuk, A. A.; Newton, M. D.; Ro¨sch, N. J. Phys. Chem. B 2003, 107, 2595. (12) LeBard, D. N.; Lilichenko, M.; Matyushov, D. V.; Berlin, Y. A.; Ratner, M. A. J. Phys. Chem. B 2003, 107, 14509. (13) King, G.; Warshel, A. J. Chem. Phys. 1990, 93, 8682. (14) Ando, K. J. Chem. Phys. 2001, 115, 5228. (15) Milischuk, A. A.; Matyushov, D. V.; Newton, M. D. Chem. Phys. 2006, 324, 172. (16) Vladimirov, E.; Ivanova, A.; Ro¨sch, N. J. Chem. Phys. 2008, 129, 194515. (17) Amber 9; Case, D. A.; Darden, T. A.; Cheatham, T. E., III; Simmerling, C. L.; Wang, J.; Duke, R. E.; Luo, R.; Merz, K. M.; Pearlman, D. A.; Crowley, M.; Walker, R. C.; Zhang, W.; Wang, B.; Hayik, S.; Roitberg, A.; Seabra, G.; Wong, K. F.; Paesani, F.; Wu, X.; Brozell, S. R.; Tsui, V.; Gohlke, H.; Yang, L.; Tan, C.; Mongan, J.; Hornak, V.; Cui, G.; Beroza, P.; Mathews, D. H.; Schafmeister, C.; Ross, W. S.; Kollman, P. A. University of California, San Francisco: 2006. (18) Wang, J.; Cieplak, P.; Kollman, P. A. J. Comput. Chem. 2000, 21, 1049. (19) Pe´rez, A.; Marcha´n, I.; Svozil, D.; Sˇponer, J.; Cheatham, T. E., III; Laughton, C. A.; Orozco, M. Biophys. J. 2007, 92, 3817. (20) Bayly, C. I.; Cieplak, P.; Cornell, W. D.; Kollman, P. A. J. Phys. Chem. 1993, 97, 10269. (21) van der Spoel, D.; Lindahl, E.; Hess, B.; Groenhof, G.; Mark, A. E.; Berendsen, H. J. C. J. Comput. Chem. 2005, 26, 1701.
5656 J. Phys. Chem. B, Vol. 113, No. 16, 2009 (22) Kutzner, C.; van der Spoel, D.; Fechner, M.; Lindahl, E.; Schmitt, U. W.; De Groot, B. L.; Grubmüller, H. J. Comput. Chem. 2007, 28, 2075. (23) Giese, B.; Amaudrut, J.; Ko¨hler, A.-K.; Spormann, M.; Wessely, S. Nature (London) 2002, 412, 318. (24) Davis, W. B.; Hess, S.; Naydenova, I.; Haselsberger, R.; Ogrodnik, A.; Newton, M. D.; Michel-Beyerle, M.-E. J. Am. Chem. Soc. 2002, 124, 2422. (25) Takada, T.; Kawai, K.; Fujitsuka, M.; Majima, T. Chem.sEur. J. 2005, 11, 3835. (26) Steinbrecher, T.; Koslowski, T.; Case, D. A. J. Phys. Chem. B 2008, 112, 16935. (27) Kurnikov, I. V.; Tong, G. S. M.; Madrid, M.; Beratan, D. N. J. Phys. Chem. B 2002, 106, 7.
Kubarˇ and Elstner (28) Voityuk, A. A. J. Chem. Phys. 2005, 122, 1. (29) Elstner, M.; Porezag, D.; Jungnickel, G.; Elsner, J.; Haugk, M.; Frauenheim, T.; Suhai, S.; Seifert, G. Phys. ReV. B 1998, 58, 7260. (30) A few of the latter specify the dependence to be λs ∝ 1- c/RDA, with the donor-acceptor distance RDA corresponding to the number of intervening nucleobases; our results are in agreement with this prediction. (31) Our result obtained with the SCC-DFTB method (ref 29). For validation, we evaluated λi for the ionization of a guanine molecule, with SCC-DFTB as well as on the PBE/TZVP and B3LYP/TZVP levels of DFT. The obtained values of 0.21, 0.20, and 0.26 eV, respectively, prove the acceptable performance of SCC-DFTB for the calculation of λi.
JP901888R
