Excited-State Properties and Relaxation Pathways of
Selenium-Substituted Guanine Nucleobase in Aqueous Solution and DNA
Duplex
The excited-state properties and relaxation mechanisms after light irradiation of 6-selenoguanine (6SeG) in water and in DNA have been investigated using a quantum mechanics/molecular mechanics (QM/MM) approach with the multistate complete active space second-order perturbation theory (MS-CASPT2) method. In both environments, the S11(nSeπ5*) and S21(πSeπ5*) states are predicted to be the spectroscopically dark and bright states, respectively. Two triplet states, T13(πSeπ5*) and T23(nSeπ5*), are found energetically below the S2 state. Extending the QM region to include the 6SeG-Cyt base pair slightly stabilizes the S2 state and destabilizes the S1, due to hydrogen-bonding interactions, but it does not affect the order of the states. The optimized minima, conical intersections, and singlet–triplet crossings are very similar in water and in DNA, so that the same general mechanism is found. Additionally, for each excited state geometry optimization in DNA, three kind of structures (“up”, “down”, and “central”) are optimized which differ from each other by the orientation of the C=Se group with respect to the surrounding guanine and thymine nucleobases. After irradiation to the S2 state, 6SeG evolves to the S2 minimum, near to a S2/S1 conical intersection that allows for internal conversion to the S1 state. Linear interpolation in internal coordinates indicate that the “central” orientation is less favorable since extra energy is needed to surmount the high barrier in order to reach the S2/S1 conical intersection. From the S1 state, 6SeG can further decay to the T13(πSeπ5*) state via intersystem crossing, where it will be trapped due to the existence of a sizable energy barrier between the T1 minimum and the T1/S0 crossing point. Although this general S2 → T1 mechanism takes place in both media, the presence of DNA induces a steeper S2 potential energy surface, that it is expected to accelerate the S2 → S1 internal conversion.
Introduction<!>System Setup<!><!>System Setup<!>QM/MM Calculations<!><!>QM/MM Calculations<!>Equilibrium Structures<!><!>Equilibrium Structures<!>Spectroscopic Properties<!><!>Spectroscopic Properties<!>Excited-State Minima<!><!>Excited-State Minima<!><!>Excited-State Minima<!>Crossing Points<!><!>Crossing Points<!><!>Crossing Points<!>Excited-State Relaxation Pathways<!><!>Excited-State Relaxation Pathways<!><!>Excited-State Relaxation Pathways<!><!>Conclusions<!><!>Author Contributions<!>
<p>In the prebiotic age, the protection of the ozone layer was not
efficient and the flux of UV radiation on Earth was much higher than
it is nowadays, resulting in an hostile atmosphere capable of causing
severe damage to DNA.1 In these conditions,
canonical purine (adenine and guanine) and pyrimidine (cytosine, thymine,
and uracil) nucleobases emerged as very photostable compounds, able
to preserve the genetic code from the deleterious effects of UV radiation.
Photostability originates from an efficient excited-state decay that
releases the excess of energy via ultrafast and efficient radiationless
deactivation processes.2−5 Today, it is well established that low-energy conical intersections
that allow the system to return to the electronic ground state in
a short time scale are responsible for the excited-state decay efficiency.6−11</p><p>Small chemical modifications can drastically alter the inherent
photostability of nucleobases.12 Particular
attention has been devoted to thiobases analogues—in which
oxygen atoms are replaced by sulfur—because, unlike their canonical
counterparts, intersystem crossing results in high quantum yields
of triplet states.13−20 For instance, 6-thioguanine (6tG) has an intersystem crossing quantum
yield of ca. 60%.17,21 The radical nature of the triplet
states render thiobases potential applications as photosensitizers
for photodynamic therapy against a variety of ailments.22 Examples are 2,4-dithiothymine,23 6-thio-2-deoxyguanine,24 or
2,6-dithiopurine.25 A natural extension
in this direction are seleno-nucleobases, where selenium replaces
oxygen or sulfur. Current studies26−29 point out that Se-nucleobases
can form stable RNA,30 DNA duplex,31 and G-quadruplex structures.32 They exhibit advantageous red-shift absorption spectra
in comparison with their thio-counterparts,33,34 and most importantly, due to their heavier atom, a faster intersystem
crossing has been observed compared to their thio-analogs—albeit
with a shorter triplet state lifetime.34</p><p>Theoretical theory is particularly well suited to unravel
the photophysical
mechanisms responsible for the efficient electronic population of
excited states.35 However, all previous
studies in sulfur- and selenium-substituted nucleobases36−43 have systematically excluded the effect of the DNA environment due
to the complexity involved, despite it could strongly influence the
photophysical properties of its chromophores.44</p><p>In this paper, we employ hybrid Quantum Mechanics/Molecular
Mechanics
(QM/MM) techniques45,46 to achieve a realistic description
of the biomolecular environment by including the effect of solvation
and that of the double helix DNA background on 6-selenoguanine (6SeG).
According to previous theoretical studies in gas phase,47 in 6SeG the S21(ππ*) state transfers its population to the S11(nπ*) state via a conical
intersection, from which the triplet states are formed. More recently,
three selenium-substituted uracils (2SeU, 4SeU, and 2,4SeU) have also
been investigated in gas phase and similar excited-state relaxation
pathways to triplet states are reported.48,49 The primary goal of this contribution is to examine whether the
main photophysical events of 6SeG are affected by the presence of
water and/or biological media.</p><!><p>The 6SeG molecule was
studied both in water solution and in a DNA environment. In both cases,
first classical molecular dynamics (MD) simulations were performed
with the Amber16 package50 in order to
obtain initial structures to use in subsequent QM/MM calculations.
Intra- and intermolecular parameters were extracted from the generalized
Amber force field (GAFF),51 while the force
field parameters of the selenium atom were taken from the parametrization
of the 2-selenouridine molecule embedded in aqueous solution.52 The RESP (Restrained Electrostatic Potential)
procedure53 was performed to get the set
of atomic charges at the B3LYP/cc-pVDZ level.54−58 In the water simulations, the TIP3P model59 was used with a system comprising one solute
molecule surrounded by 5406 water molecules inside of a truncated
octahedron simulation box of 12 Å (Figure 1a). As no constraint was applied during the
simulation, a reduced time step of 0.5 fs was employed. Initially,
the system was minimized with the steepest descent algorithm, followed
by heating in the NVT ensemble for a total time of
100 ps, scaling the temperature from 0 to 300 K. Afterward, an equilibration
regime was applied for 1 ns in the NPT ensemble in
normal conditions of temperature and pressure (P =
1 bar and T = 300 K). The final production was carried
out for 10 ns, from which the last snapshot was employed to investigate
the photophysics of 6SeG in water with a QM/MM approach.</p><!><p>QM/MM system
setup. (a) 6SeG in water solution. (b) 6SeG embedded
in DNA environment. Two different QM regions are indicated: only 6SeG
(QM1) or 6SeG paired with cytosine (QM2). All remaining atoms are
included in the MM subsystem.</p><!><p>In order to include the DNA environment, a duplex with the sequence
5′-ATGGTGCAC-3′ and 3′-TACCACGTG-5′ was
employed (Figure 1b).
This sequence was first constructed by using a nucleic acid builder
tool in Amber16 package as in previous experimental studies.60 The DNA duplex system was solvated in a water
box of 56 × 59 × 69 Å. DNA, waters, and counterions
were described with default ff99SB force fields61 and TIP3P model.59 Then, 1000
MM minimization cycles with frozen DNA were carried out, followed
by 2000 cycles without any geometric constraints. The minimized system
was heated to 300 K in a 40 ps equilibration MD simulation, followed
by 10 ns MD simulations. As in water solution, for the QM/MM study
the final MD snapshot was used, in which the oxygen atom of the fourth
G base was replaced with a selenium atom (5′-ATG-SeG-TGCAC-3′:3′-TAC-C-ACGTG-5′).</p><!><p>The QM/MM calculations
were performed using either MOLCAS8.062,63 or OpenMOLCAS64,65 interfaced with TINKER6.3.66 In aqueous
solution, the QM region contained the 6SeG molecule while all the
water molecules were described with MM. In DNA, two different QM regions
were constructed to explore effects of interbase hydrogen-bonding
interaction, see Figure 1b. The first (QM1) includes only 6SeG while the second contains the
6SeG paired with a cytosine base (6SeG-Cyt). In both models, all remaining
atoms from the DNA duplex, sugar, phosphates, water molecules, and
counterions are left in the MM region.</p><p>Regardless of the environment
(only water or solvated DNA) the QM region was described with the
state-averaged complete active space self-consistent field (SA-CASSCF)
method67 for geometry optimizations and
the multistate second-order perturbation theory (MS-CASPT2) approach68,69 on top to refine single-point energies. A total of 4 singlets and
3 triplets were included in the state-average procedure. During the
SA-CASSCF calculations, C, N, O, and H atoms were described with cc-pVDZ
and the Se atom with the aug-cc-pVDZ basis set.58,70 The active space consisted of 12 electrons distributed over the
10 orbitals; see Figure 2 for the active orbitals in water and Figures S1 and S2 for the orbitals within the DNA environment. The
MS-CASPT2 energy calculations were carried out employing larger atomic
basis sets, specifically cc-pVTZ for C, N, O and H, and aug-cc-pVTZ
for Se. In this case, the active space was augmented to 14 electrons
in 12 orbitals with σ and σ* orbitals located on the C=Se
moiety; see Figure 2.</p><!><p>Orbitals included in the active spaces used for the QM/MM calculations
in water. In green are shown the orbitals of the (12,10) active space
used in the optimizations. The two additional σ and σ*
orbitals located on the C=Se moiety are only included to refine
single point energies.</p><!><p>An imaginary level shift
of 0.2 au was employed to avoid the issues
of intruder states.71 The IPEA shift value
was set to zero.72,73 The Cholesky decomposition technique
for two-electron integrals was employed to accelerate CASSCF and MS-CASP2
calculations.74 Spin–orbit couplings
were obtained at the MS-CASPT2 level within the atomic mean-field
(AMFI) approximation.75−77 The effective spin–orbit couplings are expressed
as follows:in which ΨI and ΨJ are the perturbatively
modified electronic wave functions of the corresponding singlet and
triplet states; Hxso, Hyso, and Hzso are the x, y, and z components
of the spin–orbit operator.</p><p>We would like to emphasize
that only adiabatic states are considered
throughout this work. In this state representation, population transfers
among states of the same multiplicity are induced by nonadiabatic
couplings, which become extremely large close to conical intersections.
Population transfer between different multiplicities is mediated by
spin–orbit couplings, whose magnitude is independent of the
energy gaps. Hence, in addition to the location of relevant minimum-energy
singlet–triplet crossings, we also computed the spin–orbit
coupling matrix elements at these crossings.</p><p>For comparison,
in addition to the QM/MM calculations of 6SeG in
the presence of explicit water molecules, we also carried out computations
with the polarizable continuum model (PCM) that considers solvent
effects only implicitly.78,79</p><p>To save computational
effort, only the last saved snapshot from
the MD simulations was used as a starting point for the optimization
of excited state geometries with QM/MM. However, in order to verify
whether this approach suffices to predict deactivation pathways, we
also performed single point vertical excitation energy calculations
considering ten randomly selected snapshots taken from the classical
simulation in both environments (see Tables S1–S3). The superposition of the selected snapshots shows that the nearest
water molecules are well distributed around the 6SeG (see Figure S3), which means that the selected snapshots
are representative and give a good representation of the solvent configurations.
Since the calculated root-mean-square deviation is small for the two
lowest singlet states, with the largest deviation of 0.07 eV—within
the accepted error of the method—we conclude that the use of
one snapshot for the following discussion of the relaxation pathways
is reasonable.</p><!><p>Figure 3 shows the optimized geometries
of 6SeG in the electronic ground state (S0) in water solution
(a) and in DNA (b), as obtained with the QM(CASSCF)/MM approach. In
aqueous solution, 6SeG has a structure nearly planar with the two
hydrogen atoms of the amino group out of the six-membered ring plane
by ca. 26°. The C=Se bond length is computed to be 1.790
Å, slightly shorter than that predicted in gas phase,47 1.821 Å.</p><!><p>QM(CASSCF)/MM optimized structures of
6SeG in the electronic ground
state, in aqueous solution (a), in DNA (b), and in a pair with cytosine
(c). Selected bond lengths are in Å, and atomic numbering is
given.</p><!><p>The S0 structure of
6SeG embedded in DNA (Figure 3b) and that within the explicit
6SeG-Cyt base pair (Figure 3c) are also essentially planar. Noticeable is that the C=Se
bond lengths are 1.827 and 1.848 Å, respectively, much longer
than the corresponding C=S bond length observed in thio-counterpart
(∼1.650 Å).16,38,80 As expected, the explicit consideration of three hydrogen bonds
in 6SeG-Cyt influence the geometric parameters of the 6SeG moiety.
Accordingly, the hydrogen-bonding interaction related to the Se and
N10 atoms is much weaker than that related to the N11 atom, as suggested
by their corresponding bond lengths. The two hydrogen bonds, Se···H20
and N19···H21, are 2.635 and 2.357 Å, whereas
the O12···H22 bond length is 1.933 Å. Since the
former two hydrogen bonds are much weaker, there are small effects
on the geometric parameters around the C2–Se1 moiety. Likewise,
the stronger O12···H3 hydrogen-bonding interaction
affects the nearby bond lengths; for example, the C9–N11 and
N8–C9 bond lengths are 1.365 and 1.288 Å in 6SeG, but
1.312 and 1.336 Å in 6SeG-Cyt.</p><p>Overall, the S0 structures of 6SeG in both water and
DNA environments are similar except for the C9–N11 and C=Se
bond lengths–likely due to the pairing cytosine.</p><!><p>Table 1 collects the vertical
excitation energies of 6SeG computed in water solution and within
DNA, without and with the explicit consideration of the paired cytosine
(QM1 vs QM2). We shall first discuss the results in water solution.
The first singlet excited state, S11(nSeπ5*), corresponds to a dark state predicted at 3.13 eV in implicit solution
and 2.62 eV with explicit solvent. Relative to the S0 state,
this state is due to a single electronic excitation from the nSe nonbonding orbital localized on the Se atom, to the π5* antibonding orbital.
As mentioned above, there is a difference of around 0.5 eV between
the results obtained with the two solvation models. This difference
is expected since the degrees of freedom of the solvent are neglected
in the implicit solvation model. As known in the literature,82 lone pairs states are more affected by the explicit
solvation, especially in polar solvent, due to the formation of solute–solvent
hydrogen bonds.</p><!><p>Structure optimized with MS-CASPT2(12,10).</p><p>Structure in water optimized with
QM(CASSCF(12,10))/MM.</p><p>Structures
in DNA optimized with
QM(CASSCF(12,10))/MM: 6SeG (6SeG-Cyt).</p><!><p>The S21(πSeπ5*)
state is computed
at 3.29 eV in implicit water and at 3.51 eV with explicit solvent.
This is a singly excited configuration from the πSe orbital localized on the Se atom to the π5* orbital, and it is the bright
state (f = 0.577). It is encouraging that our calculation
with QM(MS-CASPT2)/MM level is in excellent agreement with experimental
data that reported an absorption maximum of 6SeG in aqueous buffer
solution at 357 nm (3.47 eV).34 Thus, this
electronic transition is responsible of the first peak of the electronic
absorption spectrum of 6SeG in water (see Figure S4).</p><p>The third singlet excited state, S31(πSeπ6*), is around 0.9 eV higher in energy above
the S21(πSeπ5*) state at the
QM(MS-CASPT2)/MM level of theory.
The corresponding oscillator strength is 0.016—about 1 order
of magnitude lower than that observed for the electronic transition
to the S21(πSeπ5*) state, meaning
that its contribution to the first absorption band will be small.
Accordingly, it will not be considered hereafter.</p><p>We have also
calculated the two lowest-lying triplet excited states.
They have the same electronic transition characters as the singlet
counterparts S11(nSeπ5*) and S21(πSeπ5*), i.e., T13(πSeπ5*) and T23(nSeπ5*). Their vertical
excitation energies at the QM(MS-CASPT2)/MM level of theory are 2.55
and 2.59 eV, respectively, and are overestimated with the PCM model,
particularly the T2. The next triplet state, which is described
by a T33(πSeπ6*) transition, is
situated around 1 eV above the S21(πSeπ5*) state; thus, they are not expected to be accessed from the bright
state and are not considered further.</p><p>We now proceed to discuss
the excitation energies of 6SeG in DNA.
As can be seen, both QM1 and QM2 regions (i.e., 6SeG and 6SeG-Cyt)
predict the same energetic order of the excited singlet and triplet
states and wave function character as in solution. The S1 state is spectroscopically dark with nπ* character, and the
S2 state is bright due to the ππ* character. In the S11(nSeπ5*) state the π5* orbital is delocalized
all over including the C=Se double bond. In the S21(πSeπ5*) state, the π5* orbital is the same as that in
the S11(nSeπ5*) state, while the πSe orbital is mainly composed of the p_z orbital of Se with
some minor contribution from the C3–C7 bond in either 6SeG
or 6SeG-Cyt in DNA. It is worth mentioning that the lowest-energy
states of 6SeG-Cyt are basically described by local excited states
on 6SeGua and Cyt, respectively, due to the high energy difference
between the excited states computed for each of them separately, which
are below about 3 eV for 6SeG and 5 eV for Cyt. This is confirmed
by QM(MS-CASPT2)/MM calculated vertical excitation energies on 6SeG-Cyt:
the S3 (4.64 eV) and S4 (4.77 eV) states are
local excitations on cytosine and charge transfer states, respectively,
which are on average 1.4 eV higher than S2 (see Table S6). Therefore, charge transfer states
are placed at higher energetic regions and the state character of
the 6SeG-Cyt S1 and S2 lowest states are preserved
in relation to those observed for the QM1 scheme.</p><p>Curiously,
the hydrogen-bonding interaction that appears when the
cytosine is treated quantum mechanically does not change the characters
of the S1 and S2 states as aforementioned but
it considerably affects their energies. The QM1(MS-CASPT2)/MM vertical
excitation energies of S1 and S2 states are
predicted to be 2.56 and 3.40 eV and are increased and decreased to
2.79 and 3.29 eV, respectively, with QM2(MS-CASPT2)/MM. This means
that the hydrogen-bonding interaction stabilizes the S21(πSeπ5*) state by 0.11 eV but destabilizes the
S11(nSeπ5*) state by 0.23 eV. A similar
situation is found in T13(πSeπ5*)
and T23(nSeπ5*) states. Their vertical excitation
energies are predicted as 2.35 and 2.51 eV with QM1 but 2.48 and 2.80
eV with QM2, respectively.</p><p>It is apparent that the QM(MS-CASPT2)/MM
approach predicts similar
vertical excitation energies regardless of water or DNA environment,
except for the fact that the hydrogen-bonding interactions destabilizes
the S1 and T2 states by 0.23 and 0.29 eV at
the QM(MS-CASPT2)/MM level.</p><!><p>Figure 4 shows the
lowest-lying minimum-energy
structures optimized in aqueous solution using the QM(CASSCF)/MM method.
The S1 minimum (S1-MIN) is located 2.44 eV adiabatically
above the S0 minimum at the QM1(MS-CASPT2)/MM level of
theory. It is characterized by a pyramidalization of the Se atom (by
24°) and an elongated C=Se bond length (by 0.12 Å)
with respect to the S0 minimum. The S2-MIN is predicted
at 3.34 eV above the S0 minimum. It retains the Se atom
in the molecular plane and also shows a stretched C=Se bond
(by 0.118 Å) in comparison to the S0 minimum, while
the amino group is rotated making its two hydrogen atoms almost perpendicular
to the molecular plane.</p><!><p>Minimum-energy structures of 6SeG in water solution
optimized at
the QM(CASSCF)/MM level of theory. Selected bond lengths are in Å.</p><!><p>The T1-MIN and T2-MIN also have the Se atom pyramidalized
by about
30°. The C=Se bond length is longer in the T2 minimum than in the T1 (1.941 vs 1.936 Å). These
minima are calculated 2.38 and 2.45 eV above the S0 minimum.
Experimentally, the triplet state lifetime in 6SeGua is reported to
be 835 times shorter than that observed for 6tG in aqueous solution.34 Earlier theoretical investigations in 6tG pointed
out the existence of two minimum structures in the T1 state,81 but only a single minimum is predicted for 6SeGua.
Hence, we suggest that the absence of the second minimum on T1 PES could be one of the reasons for the shorter triplet state
population.</p><p>For completeness, these minima have also been optimized
with the
PCM model. In general, the results are very similar to those obtained
with QM/MM but the pyramidalization on the Se atom becomes larger
(S1-MIN, 35°; S2-MIN, 44°; T1-MIN, 42°; T2-MIN, 41°).
The observed larger pyramidalization angles in relation to those computed
with the QM/MM method can be explained by considering how the solute–solvent
hydrogen bonds are modeled in each approach. In addition, the adiabatic
excitation energies are also influenced, with the PCM model predicting
values of 2.99, 3.18, 2.40, and 2.71 eV for the S1-MIN, S2-MIN, T1-MIN,
and T2-MIN structures, respectively. Overall, in the particular case
of 6SeG, both solvation models well agree about the relative energy
order of the minima.</p><p>The excited-state minimum-energy structures
have also been optimized
in DNA at the QM1 and QM2(CASSCF)/MM levels of theory; see Figure 5. The corresponding
S1, S2, T1, and T2 minima
can be classified into three groups, differing by the orientation
of the C=Se group–above ("up"), below
("down") or within the molecular plane ("central"),
denoted as U, D, and C, respectively. The Se atom leaves the molecular
plane pointing toward the guanine in the "up" structure,
and pointing toward the thymine in the "down" structure
(see Figure S9). As they are energetically
quasi-degenerated (Tables S3 and S4), only
the "up" conformations (referred as S1-MIN-U, S2-MIN-U,
T1-MIN-U, and T2-MIN-U) will be discussed.</p><!><p>Minimum-energy structures
of 6SeG in the S2, S1, T2, and T1 states in DNA optimized at QM1(CASSCF)/MM
(top) and QM2(CASSCF)/MM (down) levels of theory. In these structures,
the C=Se group is twisted out of the molecular plane, which
is referred to as "up" configuration and indicated by
"U" in the structure labels (see "central"
and "down" in Figures S5–S8). Selected bond lengths are in Å.</p><!><p>Conspicuously, the QM1(CASSCF)/MM calculations predict a C=Se
bond length elongated in all the structures, by more than 0.30 Å
in S2-MIN-U and by about 0.15 Å in the other minima, with respect
to the S0 minimum. In comparison, the other bond lengths
change less than 0.03 Å. The calculations with QM2(CASSCF)/MM
are consistent, with changes in the C=Se bond ranging from
0.15 to 0.32 Å and the other bond differing by less than 0.03
Å. This implies that the presence of the pairing cytosine hardly
influences the geometries of the 6SeG moiety. Interestingly, the electronic
excitations obtained with both QM regions are also very similar, as
it has been also observed in other nonnatural base pairs in DNA.83 The adiabatic excitation energies of the S1, S2, T1, and T2 are to be,
2.46 (2.58), 3.08 (2.86), 2.27 (2.29), and 2.36 (2.53) eV at the QM1(MS-CASPT2)/MM
level (QM2(MS-CASPT2)/MM in parentheses). Structures and energies
of the other conformations can be found in Figures S5–S8 and Tables S4 and S5.</p><!><p>Crossing points among
different PES are key to disentangling excited-state decay pathways.
Unlike excited-state minima, the optimization of crossing points is
much more difficult. Here, we have first optimized two-state crossing
points at the QM(CASSCF)/MM level of theory and after analyzing the
nature of the involved electronic states, a single-point energy calculation
was done with QM(MS-CASPT2)/MM, verifying whether the two states keep
degenerated. In some cases, some two-state minimum-energy crossing
points resulted in three-state degenerate points.</p><p>In water,
three minimum-energy crossing points have been found; see Figure 6. The first one corresponds
to a crossing between the S2 and S1 states (S2/S1).
This conical intersection is very similar to the S2-MIN; as a matter
of fact, both structures are quasi-degenerate in energy (S2-MIN lies
at 3.34 eV while S2/S1 is at 3.34/3.29 eV). For this reason, we computed
numerical frequencies at the SA-CASSCF level of theory (the same level
of theory employed for optimizations) in all environments, to make
sure that it is a true minimum and not simply the lowest point on
the S1/S2 intersection seam. Similar to the
three-state crossing point found for 6SeG in gas phase,47 as well as for other canonical and modified
nucleobases either in gas phase or in solution,84−87 the search for a conical intersection
between the two low-lying T1 and T2 states,
returned a three-state degenerate crossing point on which the S1 state is also degenerated. This S1/T1/T2 three-state crossing
structure has the Se atom out of the molecular plane by ca. 28°
and a C=Se bond length of 1.949 Å.</p><!><p>Minimum crossing points
of 6SeG in solution obtained with QM1(CASSCF)/MM.
Selected bond lengths are in Å.</p><!><p>The last crossing point found in aqueous solution is the singlet/triplet
crossing T1/S0, between the S0 and T1 PES. This
structure shows a distinctive large pyramidalization at the Se atom
(ca. 61°)—larger than in any of the other optimized critical
points.</p><p>Figure 7 shows the
obtained crossing points of 6SeG in DNA; as with the local minima,
only the "up" structures are shown, unless otherwise
stated (see Figures S10 and S11 for other
conformations). The two S2/S1 conical intersections
(S2/S1-U) are very similar, except for C=Se bond and its orientation
relative to the molecular plane. As expected, the presence of the
hydrogen bond shortens the C=Se bond (from 2.263 to 2.212 Å);
however, in comparison, changes in the other bond lengths can be disregarded.
The energetics of the S2/S1-U geometries are also affected by the
consideration of the base pair (see Table S5). Accordingly, while the S2 and S1 states
are at 3.28 and 3.20 eV for QM1(MS-CASPT2)/MM, they are stabilized
to 2.90 and 2.86 eV for QM2(MS-CASPT2)/MM; see Table S2.</p><!><p>Minimum crossing points of 6SeG in DNA obtained with QM1(CASSCF)/MM
(top) and QM2(CASSCF)/MM (down) levels of theory. Selected bond lengths
are in Å.</p><!><p>Singlet/triplet T1/S0 crossing points were
also optimized for 6SeG in DNA. The T1/S0-U shows a C=Se bond
length much longer than that in the S0 minimum (2.004 (QM1)/2.015
(QM2) Å in the T1/S0-U versus 1.827 (QM1)/1.848 Å in S0-MIN).
The S0 and T1 energies at T1S0-U were calculated
to be 2.58 and 2.62 eV (QM1) and 2.43 and 2.40 eV (QM2), i.e., in
contrast to the S2/S1-U conical intersection. Here, the presence of
cytosine has a negligible effect on the energies.</p><p>In addition
to the two-state intersection structures, two different
three-state intersection structures were also identified for 6SeG
in DNA: the S2/S1/T2-C and the T2/T1/S0-U. The S2/S1/T2-C crossing
structure is characterized by a much longer C=Se bond length
(1.950 and 1.926 Å) and a Se–N11–N4–N6 dihedral
angle (180° and 175°), considering QM1 and QM2 regions,
respectively. The most striking structural feature of S2/T2/S1-C,
with respect to the other intersection structures, is that the NH2 group becomes nearly perpendicular to the molecular plane.
Moreover, the C=Se bond is essentially in the molecular plane,
unlike the more pyramidalized feature of S2/S1-U and T1/S0-U. The
energies of the S2/S1/T2-C structure are 4.11/4.07/4.10 eV at the
best QM2(MS-CASPT2)/MM level of theory, indicating that this crossing
point is not accessible.</p><p>The T2/T1/S0 crossing points differ
on the orientation of the C=Se,
and accordingly, they are labeled as T2T1S0-U and T2T1S0-D. From them,
only T2/T1/S0-U is shown in Figure 7 (see Figure S10 for the
-D conformation). As can be seen, T2/T1/S0-U presents a longer C=Se
bond length (2.625(QM1) or 2.680(QM2) Å) than in the other crossing
points. The estimated energies of the T2/T1/S0 states are respectively 3.06/2.96/3.03 eV or 3.04, 3.03,
and 2.97 eV with the smaller QM1 region, showing that the hydrogen
bonding interaction has very small effects on the energetics.</p><!><p>Finally, having obtained
all the important critical points of 6SeG
in water and in DNA, we proceed now to discuss the excited-state relaxation
mechanisms, assisted by the linear interpolation internal coordinates
(LIIC) technique. Note that even when the energy barriers computed
with the LIIC technique are typically overestimated, these scans allows
us to postulate the most plausible deactivation mechanisms, which
can be further confirmed only by means of a detailed nonadiabatic
molecular dynamics study. Additionally, linear interpolation in internal
coordinates scans allow us to verify the consistency of the active
space employed for the computation of all points, whether or not there
are energy barriers due to intermediate geometries connecting the
initial and final points, and to enable a graphical presentation which
resembles the familiar potential energy surface concept. Figure 8 displays the obtained
pathways for 6SeG in water.</p><!><p>Plausible relaxation pathway of 6SeG in water
calculated at the
QM(MS-CASPT2)/MM level of theory. Relative energies of minima and
crossing points are in eV. Spin–orbit couplings at the singlet/triplet
crossings are given in cm–1.</p><!><p>According to our results, the primary UV absorption event populates
the bright S2 (1πSeπ5*) state. After
that, the system evolves barrierless along the S2 PES toward
the S2/S1 conical intersection, by dissipating 0.17 eV excess of energy.
Such an energetic profile is also confirmed using the minimum-energy
path technique with QM(CASSCF)/MM (see Figure S12). The S2/S1 conical intersection coincides with the S2-MIN
that simultaneously is also a three-state intersection, where the
S1, S2, and T2 electronic states
are quasi-degenerated. Therefore, in the vicinity of S2-MIN, the S2 electronic population can be transferred to both S1 and T2 states, through a conical intersection or a singlet/triplet
crossing, respectively. The latter process benefits from a large spin–orbit
coupling of ca. 440 cm–1 between S2 and
T2. However, which of these pathways is preferred should
be elucidated with the help of dynamical simulations.</p><p>Presuming
a partial population transfer from the S2 to
S1 state, the system can then reach the S1 minimum
(S1-MIN), located 0.90 eV below the S2/S1 conical intersection region,
as suggested by the downhill LIIC pathway. Starting from the S1-MIN
region, we investigate whether the S1/T2/T1 region, located only 0.17
eV above the S1-MIN region, could be accessed. The computed LIIC profile
does not show any barrier along this path and taking into account
that the SOC is also substantial (S1/T1 = 380
cm–1 and S1/T2 = 155 cm–1), we hypothesize that population transfer to the
T1 state should be favorable. In principle, internal conversion
to the ground state and the ISC process could occur in parallel. However,
the optimized minimum energy crossing point between the S1 and the ground state was computed to be located 1 eV above the S1
min region. Based on that, the direct relaxation pathway from the
S1 to the S0 state is suppressed, in agreement
with what is observed in other similar systems.16,48 Therefore, direct relaxation pathways to the S0 were
not considered by us.</p><p>The effects of sulfur-to-selenium substitution
in the photophysics
promotes a substantial decrease of the triplet lifetime, as observed
by the experiment.34 Our LIIC scans shows
that the ISC process is favorable for 6SeG in water; however, from
static calculations, we are not able to quantify how much more efficient
this process is in relation to 6tG. Thus, only nonadiabatic dynamics
simulations could be used in order to better understand why the 6SeG
triplet lifetime is shorter. From here, there is a barrierless pathway
toward the T1-MIN (0.26 eV below the three-state S1/T2/T1 crossing
region), from which radiative and nonradiative processes can occur.
The system can emit phosphorescence at 1.63 eV. Alternatively, as
shown by the LIIC pathway, the singlet/triplet T1/S0 crossing located
adiabatically 0.26 eV above the T1 minimum, allows an intersystem
crossing to the S0 state. The calculated spin–orbit
coupling for this crossing (283 cm–1), despite smaller
than the previous values, could be sufficient for an efficient deactivation.
Although, with MS-CASPT2 the energy gap at the T1/S0 crossing is still
0.33 eV at the CASSCF optimized structure, we assume that population
to the ground state is viable. For comparison, in gas phase, the T1/S0
crossing point is predicted 0.11 eV above the T1-MIN at the MS-CASPT2
level of theory47 – a little lower
than that in aqueous solution.</p><p>The excited-state deactivation
mechanisms of 6SeG (QM1) and 6SeG-Cyt
(QM2) in DNA have been also investigated, see Figure 9. In DNA, three possible deactivation mechanisms
can be followed from the Franck–Condon region. One of them
involves a barrierless relaxation to the S2 minimum with
the up configuration i.e. S2-MIN-U, by releasing
0.32 eV of energy, as computed at the QM1(MS-CASPT2)/MM level (Figure 9a). Once in the S2-MIN-U
region, the system evolves toward an internal conversion region with
the S1 state (S2/S1-U, with the up configuration
too, which is located adiabatically 0.20 eV above the S2-MIN-U region.
From here, population is transferred to the S1 state, finally
reaching barrierlessly the S1 minimum S1-MIN-U, in which
the up configuration is kept. Analogously, compared
to the 6SeG in water, the S1 of 6SeG and 6SeG-Cyt in DNA
need to overcome 1.03 and 0.67 eV energy barriers to reach the S1/S0
conical intersection region, whereas LIIC paths connecting the S1 and triplet manifolds are nearly barrierless; consequently,
we believe that S1 directly deactivating to the ground
state is inefficient. The S1 minimum region also represents
a three-state crossing region among the S1, T2, and T1 electronic states, on which the population can
be transferred from the S1 to T1 or T2 states, via intersystem crossings (S1/T1 or S1/T2). The S1 → T1 intersystem crossing process is more likely
due to its larger S1/T1 spin–orbit coupling
of 439.4 cm–1, in accord to the El-Sayed rule.91 Once on the T1 state, the system
evolves to T1-MIN-U, from which the intersystem crossing with the
S0 state, T1/S0 (T1/S0-U, SOC = 273
cm–1), can be reached after overcoming an energetic
barrier of ca. 0.30 eV. Due to this barrier, the system will be trapped
a while in the T1 state, before hopping to the S0 state.</p><!><p>Plausible relaxation pathway of 6SeG in DNA calculated at the (a)
QM1(MS-CASPT2)/MM and (b) QM2(MS-CASPT2)/MM levels of theory. Relative
energies of minima and crossing points are in eV. Spin–orbit
couplings at the singlet/triplet crossings are given in cm–1.</p><!><p>Alternatively, the S2 state can also relax toward another
S2 minimum, the S2-MIN-D (see Figure S13). In such a case, and as discussed previously, after reaching
the T1-MIN-D, it is necessary to overcome an energy barrier of 0.26
eV to access the T1/S0-D crossing point. Finally, if instead, the
S2 state relaxes to the S2-MIN-C minimum, its further relaxation
to the S1 state via the S2/S1/T2-C three-state intersection
point encounters a comparable barrier in the S2 state because
this process involves a large rotation of the N11H2 group
(see Figure S14). Accordingly, this pathway
does not compete with the other two. However, if S2-MIN-C can be first
converted either into S2-MIN-U or S2-MIN-D in the S2 state,
the system initially populated at S2-MIN-C still can decay to the
T1 state efficiently.</p><p>The calculations of the excited-state
relaxation pathways of 6SeG
in DNA using the GC base pair (Figure 9b) are similar to those described above, with few exceptions.
The S2 potential energy is much steeper, so that the initial
S2 relaxation pathway toward its minimum, S2-MIN-U, should
be faster. The internal conversion process from the S2 to
S1 states via S2/S1-U is essentially barrierless (0.02
eV here versus 0.20 eV without C). This is mainly attributed to the
hydrogen bond interaction between the 6SeG and C. Intermolecular hydrogen
bonds stabilize the 6SeG ππ* excited
state level, which significantly reduces the excited state lifetime.
In fact, several similar phenomena were reported in the literature.88−90 For example, Temps et al. observed that the formation of guanine
and cytosine Watson–Crick base pairs results in ultrafast nonradiative
transition.88 Around the S1-MIN-U minimum,
the system will undergo intersystem crossing to both T2 and T1 states. But, as mentioned above, the S1 → T1 intersystem crossing is preferred because
it has a larger S1/T1 spin–orbit coupling
of 399 cm–1, justified by the El-Sayed rule.91 Once in the T1 state, the barrier
separating T1-MIN-U from T1/S0-U increases to 0.14 eV, i.e., 0.14
eV higher than that without C. For completeness, we mention that when
the S2 system relaxes to its other minimum S2-MIN-D, similar
excited-state relaxation pathways are identified (see Figure S15) and the corresponding pathways from
S2-MIN-C are not studied here because of the high energies of the
S2/T2/S1–C triplet crossing point, which makes it unlikely
(see above).</p><p>Figure 10 displays
a schematic summary of our proposed relaxation mechanisms. In short,
the environment surrounding the 6SeG greatly affects how fast the
internal conversion from the bright state to the lowest singlet state
takes place. A high energy barrier is computed in the gas phase,47 disappearing in water. Depending on which structure will
be primarily accessed in DNA, the population transfer to the S1 PES may be faster or slower. Afterward, the population will
be easily transferred to the triplet states via ISC, since the internal
conversion to the S0 is less important.</p><!><p>Possible deactivation
mechanisms of 6SeG in gas phase, water, and
DNA.</p><!><p>A multiscale QM/MM strategy including the very accurate MS-CASPT2
method for the QM part has been employed to explore the excited-state
properties and excited-state decay pathways of 6SeG in water and in
a DNA environment, and compare these to those in gas phase. In the
DNA environment, the effect of the pairing cytosine was also investigated
by including it into the QM region. Spectroscopically, in all the
three environments, 6SeG is characterized by the same two lowest excited
singlet states, i.e. S11(nSeπ5*) and S21(πSeπ5*), corresponding to a dark and the spectroscopically
bright excited states, respectively. Two triplet electronic states,
i.e., T13(πSeπ5*) and T23(nSeπ5*), were found to be energetically lower than
the bright S2 state. We find that the order of the singlet
and triplet states and their characters does not depend on the environment;
moreover, the vertical excitation energies in water and DNA are very
similar to each other, as long as only the 6SeG is considered quantum
mechanically. The inclusion of explicit hydrogen-bonding interactions
stabilizes the S21(πSeπ5*) state and destabilizes
the S11(nSeπ5*) state.</p><p>For comparison,
the inclusion of solvent is also taken into account
by implicit methods to find that PCM model overestimates both S11(nSeπ5*) and T23(nSeπ5*) states, in relation to the values computed with the QM/MM method.</p><p>In water, we have found that the excited-state relaxation pathway
corresponds to an efficient relaxation from the initially populated
S21(πSeπ5*) state to the
T13(πSeπ5*) state. After
irradiation, the S2 state quickly relaxes from the Franck–Condon
region to an S2 minimum, near which a S2/S1/T2 three-state intersection region induces both
S2 →S1 internal conversion and S2 → T2 intersystem crossing processes. Subsequently,
it is expected that both the S1 or T2 states
decay further to the T1 state via either S1 →
T1 intersystem crossing or T2 → T1 internal conversion. Both intersystem crossing processes
benefit from substantial spin–orbit couplings (S1/T2 = 440 cm–1 and S1/T1 = 380 cm–1). Once
in the T1 state, the system should be trapped for some
time, due to the existence of a small energy barrier between the T1 minimum and the T1/S0 crossing point,
estimated to be 0.26 eV.</p><p>The inclusion of DNA environment increases
the complexity of the
excited-state relaxation pathways of 6SeG due to the existence of
different conformational structures, that depend on whether the C=Se
group is oriented above ("up"), below ("down")
or within ("central") the molecular plane. In general,
it is possible to identify a preferred pathway, with resemblances
to that in water. For example, starting from the up configuration,
the excited S2 system first arrives at its S2 minimum, and then, a nearby S2/S1 conical
intersection drives the system down to the S1 state. An
S1 → T1 intersystem crossing process
then occurs in the vicinity of the S1 minimum, and an energy
barrier of 0.30 eV traps the system in the T1 state before
it can go to the S0 state. A similar excited-state decay
pathway is found from the down configuration, while the central configuration
is predicted not to be efficient because of a large energy barrier
in the S2 state. Which of these pathways is most probable
and in which amounts could be predicted by dynamical simulations.</p><p>The inclusion of the pairing cytosine has little influence on the
obtained relaxation mechanism, except that the S2 →
S1 internal conversion processes is predicted to be faster
and larger energy gaps between the T1 minima and the T1/S0 crossing points are found.</p><p>The present
work provides an important step toward understanding
the excited state properties of chemically modified nucleobases and
its response to light, on the way to find better photosensitizers,
for example for phototherapies.</p><!><p>Minima and intersection
structures for the other structures
optimized in DNA, excited state relaxation paths for the down (D)
and central (C) structures optimized in DNA, simulated electronic
absorption spectrum in water, minimum energy path in water, and Cartesian
coordinates of the optimized structures (PDF)</p><p>jp0c10855_si_001.pdf</p><!><p># Y.-G.F. and D.V. contributed equally to this work.</p><!><p>The authors declare no
competing financial interest.</p>