Proton-Transfer Reactions in Reaction Center of Photosynthetic

Jun 4, 2009 - Department of Chemistry, Graduate School of Science, Nagoya University, Furocho, Chikusaku, Nagoya 464-8602, Japan, Department of ...
0 downloads 0 Views 5MB Size
J. Phys. Chem. B 2009, 113, 8993–9003

8993

Proton-Transfer Reactions in Reaction Center of Photosynthetic Bacteria Rhodobacter sphaeroides Yu Kaneko,†,⊥ Shigehiko Hayashi,‡,¶ and Iwao Ohmine*,†,§ Department of Chemistry, Graduate School of Science, Nagoya UniVersity, Furocho, Chikusaku, Nagoya 464-8602, Japan, Department of Chemistry, Graduate School of Science, Kyoto UniVersity, Kitashirakawa, Sakyo-ku, Kyoto 606-8502, Japan, and Fukui Institute for Fundamental Chemistry, Kyoto UniVersity, Nishihiraku-machi 34-4, Sakyo-ku, Kyoto 606-8103, Japan ReceiVed: January 30, 2009; ReVised Manuscript ReceiVed: April 30, 2009

The mechanism of proton uptakes by the secondary ubiquinone (QB) in the reaction center of the photosynthetic bacteria Rhodobacter sphaeroides is investigated theoretically. Two protons are transferred to the secondary ubiquinone (QB) upon two electron transfers to it. The pathways of the proton transfers are explored through molecular dynamics and free-energy calculations by a molecular mechanical method and potential energy surface calculations by a hybrid ab initio quantum mechanical/molecular mechanical method. Initial locations of donor protons are identified at Glu-L212 and at Asp-L210. It is shown that the first proton transfer takes place from Glu-L212 to QB through their direct hydrogen bond only when the second electron is transferred to QB. The second proton transfer from Asp-L210 to QB proceeds through long-range hydrogen-bond network bridges connecting them. The hydrogen-bond network bridges are occasionally formed with several water molecules in the water molecular structural fluctuation dynamics of the “proton channel”. The activation energy barrier along the second proton-transfer path thus formed is consistent with the experimental rate. It is also found that there exist very strong interactions among water molecules and a protonated carboxyl group of Asp-L210, suggesting formation of a hydronium ion in the surroundings of negatively charged acidic groups in the proton channel. 1. Introduction The bacterial photosynthetic reaction center of Rhodobacter sphaeroides functions as an enzyme for reduction of an ubiquinone10 with utilization of light energy in bacterial photosynthetic activities.1-6 The reaction center (RC) is a membrane-integrated complex which consists of three transmembrane subunits (chain L, M, and H), four bacteriochlorophylls, two bacteriopheophytines, one nonhem iron, and two ubiquinones, called the primary quinone, QA, and the secondary quinone, QB. Figure 1a summarizes chemical events in the photocycle. Two electrons are transferred through QA to QB upon absorptions of two photons, and then, QB accepts two protons from cytoplasmic solution. After this full reduction, QB in a quinol form is released from the protein. There have been numerous experimental and theoretical studies on the photochemical process of the important reactions in biology. The electron-transfer mechanism of the RC unit up to QB has been studied experimentally by flash photolysis methods, and the electronic states and their couplings involved in the electron transfer up to QA have been investigated theoretically by QM and QM/MM calculations.7-10 Significant amounts of knowledge are being accumulated for this electronic part of the reaction. * To whom correspondence should be addressed. E-mail: ohmine@ aqua.chem.nagoya-u.ac.jp. Phone: +81-52-789-2481. Fax: +81-52-7893551. † Nagoya University. ‡ Graduate School of Science, Kyoto University. § Fukui Institute for Fundamental Chemistry, Kyoto University. ⊥ E-mail: [email protected]. ¶ E-mail: [email protected].

On the other hand, a detailed mechanism of two proton uptakes by QB, which occur after the two electron transfers, has not been fully explored yet. For example, how protons merge into the protein from cytoplasmic solution and at which part of the protein they are located in various QB states are not wellknown. The entrance area of two incoming protons has been considered to be around Asp-H124, His-H126, and His-H128.11-15 From this area, as various mutagenesis experiments indicate, the protons migrate to the inside of the protein through a negatively charged “proton channel” consisting of Asp-L210, Asp-L213, Ser-L223, Asp-M17, and Glu-H173.1,5,11,12,16-26 Figure 1b depicts the molecular structure of the proton channel. The incoming protons likely locate temporarily at some of these acidic residues. Several FTIR analyses have been made to specify the protonation positions among the acidic groups and Glu-L212 in various QB states. Glu-L212 has been shown to be protonated in the first electron-transferred state (QAQB-). Some FTIR experiments reported a transient protonation of Asp-L21027,28 in addition to the protonation of Glu-L212 in this state, while other FTIR experiments showed protonation of none of those acid residues.6,29-31 It is, however, natural to assume the second excess proton exists in the proton channel even in the QAQB- state. Otherwise, the long-range proton migration from bulk solution to some amino acid residues in the proton channel has to take place only after the second electron attaches to QB. Such a migration would hardly occur on a time scale experimentally observed. The position and dynamics of this second proton will be discussed in the present report. When the second electron attaches to QB (QAQB2-), two protons are transferred to QB through the proton channel, as the mutagenesis experiments have shown. The exact mechanism of these proton transfers has not been fully explored yet due to

10.1021/jp9008898 CCC: $40.75  2009 American Chemical Society Published on Web 06/04/2009

8994

J. Phys. Chem. B, Vol. 113, No. 26, 2009

Kaneko et al.

Figure 1. Scheme of the electron and proton transfers in the reaction center. (a) The overall scheme. (b) The proton-transfer pathways to QB in the proton channel. Blue and red arrows represent the electron-transfer and the proton-transfer pathways, respectively. BCl and BPh are bacteriochlorophyll and bacteriopheophytine, respectively.

experimental difficulties to observe the transient protonation; for example, a direct monitoring of stable or transient protonation of these acidic groups by FTIR analysis has not been performed at the second electron-transfer process. In the present study, we perform theoretical analyses on the mechanism of the QB proton uptake after the second electron transfer by using theoretical methods, that is, molecular dynamics (MD) simulations and free-energy analysis using a molecular mechanics (MM) method and the potential energy surface calculation using a hybrid ab initio quantum mechanical/ molecular mechanical (QM/MM) method. Locations of the excess protons on the acid residues and the reaction paths of the proton transfers from the protonated residues to QB are investigated. Section 2 describes the theoretical methods used in the present study. Results and Discussion are given in section 3. Conclusions are given in section 4. 2. Methods Section The initial protein system is constructed based on a crystallographic structure (PDB code, 1AIG).3 The Amber parm99 force field parameters32 are used for the protein. To describe the water solvation around the proton channel, 238 water molecules are added to the system around the Nδ1 atom of HisH126 with a 17 Å cap potential. A single water sphere of five to six water layers is applied to treat water molecules around the entrance area. The X-ray crystallographic structure (its resolution is ∼2.0 Å) is used to provide the initial locations of crystal water molecules inside of the protein. For a treatment of ambulant water molecules in the proton channel and the surface area, we performed long equilibrium MD simulations (∼20 ns) by applying a cap water potential. The SPC/E model is employed for the water model.33 The total number of atoms in the system is 14991. In order to examine the proton transfer to QB- after the second electron transfer, QA-QB-, QAQB2-, QAQBH-, and QAQBH2 states are considered. Atomic partial charges of the ring moiety of ubiquinone10 at the states Q, Q-, Q2-, QH-, and QH2 are determined with ESP-derived charges in isolated condition at the HF/6-31G+(d) level. Atomic charges of bacteriochlorophyll,

bacteriopheophytine, and the tail part of ubiquinone10 are determined with ESP charges at the HF/6-21G level. The Amber parm99 force field parameters32 are used for the internal and Lennard-Jones potential functions. MD Simulation using MM Potentials (MM/MD). MD equilibrium simulations for the protein with the cofactors are performed by Amber version 6. Since the proton-transfer process is considered to be localized in the regions of QB and the proton channel, a partial protein region including groups within 17-20 Å of the proton channel is only allowed to move (see Figure S1 in Supporting Information). No periodic boundary condition is used, and a 15 Å cutoff is applied for nonbonding interactions. Temperature is controlled to be 300 K by the Berendsen algorithm.34 Bonds including hydrogen atoms are fixed by the SHAKE algorithm.35 The time step used is 1.0 fs. Thermodynamic Integration Method. Free-energy differences between the protonated states of Asp-L210 and Asp-L213, that is, [Asp-L210-H, Asp-L213-] and [Asp-L210-, AspL213-H], are evaluated by Gibbs program of Amber version 7. Dummy atoms are introduced to switch the protonation states of aspartates. Only the electrostatic potential is switched between the ionized state and the protonated state. The free-energy contributions originating from changes in the nonbonding interactions of the aspartates with the surrounding groups upon the protonation and deprotonation are evaluated by gradually switching the force field parameters with a standard thermodynamic integration procedure. The contributions of the internal potential functions are canceled by simultaneous changes of the protonation states of Asp-L210 and Asp-L213. The initial structures of QA-QB- and QAQB2- with the [Asp-L210-H, AspL213-] state are taken from structures after 3.0 ns of equilibration. The free-energy perturbation calculations are performed in forward and backward directions for 1.0 ns each. An equilibration calculation is carried out for 1.0 ns before the backward calculation. The average of the free-energy changes obtained by the forward and backward simulations is used as the value of the free-energy difference.

Reactions in Reaction Center of Rhodobacter sphaeroides QM/MM Method. The QM/MM method previously employed36 is implemented in the GAMESS program package.37 The QM and MM regions are treated by B3LYP/6-31G(d,p) and the Amber parm99 force field, respectively. Atoms contained in the QM region are defined in the Results and Discussions below. Dummy hydrogen atoms as the link atoms are set in the QM region. No cutoff is used for nonbonding interactions. Structures at 0 K obtained by simulated annealing MD simulations are used for the initial structures of the QM/ MM optimization. A SHAKE algorithm is employed for bond length constraints in “minimum-energy path” calculations. The minimum-energy path means here a path which has the lowest potential energy along a given reaction coordinate. It provides a rough overall potential profile of the proton transfer. We have attempted several reaction coordinates to determine a minimumenergy path of a long-range proton transfer. The minimumenergy pathway must smoothly connect between the initial and the final state. Definition of Hydrogen Bond (HB). Given an O1-H1 bond of a water molecule and an O2 atom of another water molecule, a HB is defined to be formed between O2 and H1 if the distance between O1 and O2 is less than or equal to 3.0 Å and the exterior angle of O1-H1-O2 is less than or equal to 30°. The protontransferable pathways temporarily formed between a proton of Asp-L210 and QBH- (see below) are identified in sampled configurations from the 20 ns equilibrium trajectory of the QAQBH- state. Pathways consisting of 10 and more intervening water molecules are discarded. 3. Results and Discussions The positions of excess protons and their transfer to the second quinone, QB, are investigated for two states, namely, QA-QB- and QAQB2-. It is now well-established that Glu-L212 is protonated after the first electron transfer.38-40 One proton transfer is thus expected to take place from this amino residue to QB when the second electron is attached to QB. As the position of another excess proton, several residues such as Asp-L210, Asp-L213, and Asp-M17 have been suggested experimentally.1,16,19-21,24,29,31 Among them, protonation of Asp-M17 is less probable because the group is located at protein surface and well solvated by bulk water molecules. AspL210 and Asp-L213 are thus the candidates for the proton donor of the QB proton uptake. We first evaluate the relative proton affinity of Asp-L210 and Asp-L213 by a thermodynamic integration method using MM/ MD simulations. Proton-transfer pathways from these protondonor groups to QB are then determined by MM/MD simulations. The potential energy profiles of these proton-transfer paths are evaluated by QM/MM calculations. Evaluation of Protonation States of Aspartate Groups in the Proton Channel. The region around Asp-L213 and AspL210 is electrostatically very negative, and one of these aspartates is expected to be protonated in the QA-QB- state. The pKa’s of the carboxyl groups of these residues can be largely altered by water penetration into the proton channel.41-43 In the MM/MD simulations for the QA-QB-, two distinct states regarding the water penetration, that is, low and high water occupancy states, are found. Figure 2a and b depicts snapshots of those states. A MM/MD simulation starting from the X-ray crystallographic structure of the QA-QB- state results in the low water occupancy state in which about three water molecules occupy the channel. The high water occupancy state is obtained by a rotation of the carboxyl group of Asp-L210 around the CR-Cβ bond induced by external forces, which enhances the

J. Phys. Chem. B, Vol. 113, No. 26, 2009 8995

Figure 2. Snapshots of amino acid and water structures around the “proton channel” in the QA-QB- state observed in the MM/MD trajectory. The structures of [Asp-L210-H, Asp-L213-] in the low and high water occupancy states are shown in (a) and (b), respectively. The structures of [Asp-L210-, Asp-L213-H] in the low and high water occupancy states are shown in (c) and (d), respectively. Water molecules in the “proton channel” are those in green circles. All of the 3D modeling snapshots of the protein structures are generated by VMD.44

rearrangement of the hydrogen-bond network and thus allows water migrations from bulk solvent in the MM/MD simulations. Once the migration of water molecules is achieved, this high water occupancy state with five to eight water molecules is stable during a 20 ns equilibrium simulation even in the absence of the external forces. Those high and low water occupancy states are also found for QAQBH-. Mutual thermodynamic stability of the protonated states among the aspartates is examined for the QA-QB- state. We compute the free-energy difference between two protonated states, that is, the state in which Asp-L210 and Asp-L213 are protonated and unprotonated, respectively, that is, [AspL210-H, Asp-L213-], and the other state [Asp-L210-, AspL213-H] with the alternated proton configuration, by using the thermodynamic integral method with a switching function in the MM/MD simulations. It is found that the former state, [AspL210-H, Asp-L213-], is more stable by 35.1 (20.5) kcal/mol than the latter one, [Asp-L210-, Asp-L213-H], in the high (low) water occupancy state. The preference of the ionized form of Asp-L213 in the QA-QB- state is consistent with experimental observations.6,29 The side chains of the aspartates form in very different conformations for the different protonated states. Figure 2c,d depicts the side chain conformational changes in the proton channel region. In the [Asp-L210-H, Asp-L213-] state, a strong salt bridge is made between the negatively charged Asp-L213and a proximal arginine, Arg-L217, contributing to the large stabilization of this protonation state. On the other hand, such a salt bridge cannot be formed in the [Asp-L210-, Asp-

8996

J. Phys. Chem. B, Vol. 113, No. 26, 2009

Kaneko et al.

L213-H] state, where Asp-L213 makes a strong HB with negatively charged Asp-L210 instead. The free-energy difference between these protonated Asp states is also computed by the same method for QAQB2-. The state of [Asp-L210-H, Asp-L213-] is again found to be more stable by 27.3 (18.5) kcal/mol than the state of [Asp-L210-, Asp-L213-H] in the high (low) water occupancy state. The large stabilization of the [Asp-L210-H, Asp-L213-] state is also due to the formation of a strong salt bridge between AspL213- and Arg-L217 in QAQB2-. We will mainly deal with the high water occupancy state hereafter because the formation of the HB network in the high occupancy state makes the proton transfer along the network possible, whereas such a HB connection cannot be formed in the low water occupancy state. The Proton Transfer from Glu-L212 to QB-. The protontransfer processes to QB are induced by the second electron transfer to QB. Glu-L212 is known to be protonated in QAQBand QA-QB- states. It is found from the MM/MD simulations that protonated Glu-L212 forms a stable direct HB with QB- in QA-QB- and QAQB2- states. The proton transfer is thus expected to take place first from Glu-L212 to QB on this HB. We perform QM/MM calculations to examine the potential energy profile along this proton-transfer path. Figure 3a depicts the QM region, which consists of the ring of QB-, the side chain of Glu-L212, and a water molecule. The structures taken from the equilibrium MM/MD simulations (see the Methods Section) for QA-QBand QAQB2- states are optimized by the QM/MM method. In QA-QB-, the optimized structure is similar to those observed in the MM/MD trajectory. The proton is initially placed on GluL212 in the QM/MM optimization, and a minimum-energy path along a reaction coordinate, C1, from the initial state is determined. The reaction coordinate C1 is defined as

|





| |





C1 ) RO1 - RH1 - RO2 - RH1

|

(1)

Figure 3b shows the computed potential energy profile along this coordinate in QA-QB-. One can see in the figure that the potential energy increases monotonically by more than 18 kcal/ mol when the proton moves from Glu-L212 to QB-. This shows that the proton transfer does not proceed in the QA-QB- state, which is consistent with experimental evidence.1,38,45-47 The overall QM/MM potential energy can be decomposed into the energies of the QM and MM regions and the interaction between them. It is seen in Figure 3b that both QM and QM/ MM interaction energies increase monotonically while the MM potential energy is almost constant along the reaction coordinate. This means that the changes in the local electronic nature of the QB- and Glu-L212 and the concomitant change of their interaction with the surrounding protein are the main factors to determine this potential energy profile. Next, we compute the minimum-energy path of the proton transfer in the QAQB2- state, that is, after the second electron transfer. One can see in Figure 3c that the potential energy profile along the minimum-energy path is strikingly different from that in the QA-QB- state (Figure 3b). The potential energy decreases monotonically by ∼13 kcal/mol along C1, indicating that the proton is smoothly transferred to QB from Glu-L212 upon the second electron transfer. As seen in Figure 3c, the energy contribution of the QM region mainly determines the overall QM/MM potential energy profile; the QM energy steeply decreases by more than 23 kcal/mol along C1 because of a large proton affinity of QB2-. On the other hand, the interaction between the QM and MM regions increases by ∼10 kcal/mol

Figure 3. Proton transfers from Glu-L212 to QB in QA-QB-, and QAQB2- states. (a) The snapshot of the defined QM region. The ring part of QB, the side chain of Glu-L212, and a water molecule are defined as the QM region. (b) The potential energy profile of the proton transfer between QB and Glu-L212 in QA-QB-. (c) The potential energy profile in QAQB2-. Blue, red, green, and cyan lines indicate the energies of the overall system, QM, QM/MM interaction, and MM, respectively. (d) Schematic potential energy profiles of the first proton transfer from Glu-L212 to QB coupled with the second electron transfer. The relative energy of the QB- and QB2- states, not estimated in the present calculation, is arbitrarily chosen.

upon the proton transfer. The drastic change of the potential energy profile of the proton-transfer process before and after the second electron transfer clearly illustrates a tight coupling between the electron- and the proton-transfer processes. For the sequential order of the second electron transfer and the first proton uptake of QB-, two kinds of mechanisms, namely, (1) a proton-activated electron-transfer mechanism and (2) a concerted proton- and electron-transfer mechanism, have been proposed. There have been several experiments, for example, a driving force assay,48 to find this sequential order, but so far, they could not have distinguished between these two mechanisms. The monotonically increasing and decreasing

Reactions in Reaction Center of Rhodobacter sphaeroides potential slopes (QB- and QB2- in Figure 3d) give rise to the electron-transfer transition coupled with the proton transfer in a concerted manner. Although the present calculation cannot directly address this issue, the above result suggests a concerted proton- and electron-transfer mechanism. There is, of course, a possibility of the protonated Glu-L212 being more strongly hydrated than the configuration presented here; then, the protonated Glu-L212 does not make the direct hydrogen bond with quinone but rather forms a hydrogen bond bridge to the latter through one or two water molecules. In such a case, the rate of the proton transfer from Glu-L212 to quinone could be slow. The Proton Transfer to QBH- from Asp-L210. After the first proton uptake of QB upon the second electron transfer, QBH- further accepts a second proton to form QBH2 in several tens of milliseconds.49 As described above, Asp-L210 is a very good candidate for the proton donor of the second proton uptake of QBH-. We investigate the mechanism of the proton transfer from Asp-L210 to QBH-. Since Asp-L210 is located at a distance from QBH- (Figure 1b), a long-range HB network needs to be formed between these proton donors and acceptors in order to make the proton transfer possible. We first examine the formation of such a hydrogen-bond network in the proton channel moiety. A 20 ns equilibrium trajectory is obtained by the MM/MD calculation (see the Methods Section), and the formation of the transient proton transferable pathways is identified among conformations sampled from the trajectory. The definition of the HB described in the Methods Section is employed. Figure 4a shows relative populations of various protontransferable pathways, consisting of 3-10 water molecules, formed in the sampled conformations. The relative population of the pathways is defined as the ratio of the number of MM/ MD-sampled configurations making the HB bridges to the number of the total configurations. Note that the relative populations are less than 0.01, showing that the protontransferable pathways are seldom formed. The pathway composed of seven water molecules has the highest probability among the transferable pathways. The pathway with five water molecules yields the second highest probability. The protonated carboxyl group of Asp-L210 is rotating around the dihedral angle (Oδ2-Cγ-Cβ-CR), where Oδ2 is the oxygen atom attached to a proton. One can see in Figure 4b that this carboxyl group yields the highest probability, being at 100-120° and the nest at 180°. The orientation of the protonated carboxyl group of Asp-L210 plays an important role in the formation of the proton-transferable pathways. Figure 4b shows that there exists the strong correlation between the population of the proton-transferable pathway and that of the orientation of the carboxyl group in the sampled conformations of the trajectory; the population of the proton-transferable pathway is large when the distribution of the dihedral angle is large, indicating that the proton-transferable pathways are formed at stable orientations of the carboxyl group. A shoulder is seen for the population of the proton-transferable pathway at 100-120°, where the population of the dihedral angle is very large. A small peak is also observed in the distribution of the protontransferable pathway at -75°, where the peak of the dihedral angle population also exists. Figure 4c identifies the correlation between the number of water molecules involved in the transient HB bridges and the carboxyl group orientation. The pathway at the peak distribution around 180° is composed of 7 water molecules. The shoulder and the peak in the distribution of the pathway around 100-120

J. Phys. Chem. B, Vol. 113, No. 26, 2009 8997

Figure 4. Relative population of the proton-transferable pathway from Asp-L210 to QBH- defined as the ratio of the number of MM/MD sampled configurations making the HB bridges to the number of the total configurations. (a) The population of the number of water molecules in the identified pathway. (b) The distribution (blue line) of the dihedral angle (Oδ2-Cγ-Cβ-CR) of Asp-L210, where Oδ2 is the oxygen atom attached to a proton, and the relative population (red line) of the pathways against the dihedral angle. The total distribution of the dihedral angle is 1.0. The total population of the pathway formed is 0.015. (c) Two-dimensional plot of the relative population of the pathways against the number of water molecules participating to make the bridge (the vertical axis) and the dihedral angle (the horizon axis).

and -75° are characterized by the pathways with five water molecules. Among them, the paths with five water molecules at 100-120° play the most important role in the second protontransfer process, as seen below. Figure 5a and b shows a HB pathway formed with seven water molecules at 180° and that with five water molecules at 100°, respectively. These two pathways have mutually the almost same HB connections up to the third water molecule from QB, that is, W1-W3. At 100°, where the system has the largest population (see Figure 4b), the protonated carboxyl group of Asp-L210 is stabilized by making a firm HB with W5, and this W5 makes stable HBs with W4. This HB network from the protonated carboxyl group is usually terminated at W4 (see Figure 6a). W4, however, undergoes a large flip on some occasions (with the probability of 0.13) to form HB with W3, and then, a HB bridge passing through W1-W5 is created (Figures 6b and 5a). At 180°, the protonated carboxyl group of

8998

J. Phys. Chem. B, Vol. 113, No. 26, 2009

Kaneko et al. below), from the R state are determined by the QM/MM calculations. The QM/MM-optimized structure of the R state is depicted in Figure 7b. It is noteworthy that the optimized O-O distance between Asp-L210 and W5 is remarkably short, 2.59 Å, indicating the existence of the very strong HB between them. A minimum-energy path from the R state along a reaction coordinate, C2

|





| |





C2 ) RO - RH2 - RO - RH1 Figure 5. Snapshots of the identified proton-transferable pathways in the QAQBH- state. (a) A structure of the pathway consisting of five water molecules. (b) A structure of the pathway consisting of seven water molecules.

Figure 6. Snapshots of proton-transferable pathway formation by the hydrogen bond network interchange inside of the proton channel. (a, b) A pathway with five water molecules. A pathway with five water molecules is bridged between the bond between W3 and W4 when W4 flips to the right torsional angle (see text), (a) before the W4 flipping and (b) after the W4 flipping. (c, d) A pathway with seven water molecules. A proton-transferable pathway with seven water molecules is formed when the hydrogen bond between Asp-L210 and W7′ is established by the rotation of the protonated carboxyl group of AspL210, (c) before the rotation and (d) after the rotation. Panel (b) corresponds to Figure 5a, and panel (d) corresponds to Figure 5b.

Asp-L210 usually making a HB with W6 occasionally (with the probability of 0.07) rotates around the dihedral angle, H-Oδ2-Cγ-Cβ, to make the new HB bond with W7′ (Figure 6 c,d), and then, a HB path involving seven water molecules, W7′, W6′, W5′, W4′, and W3-W1, can be established (Figures 6 and 5b). We then examine the potential energy profiles of the proton transfers along the pathways consisting of five and seven water molecules by the QM/MM calculations. The QM/MM optimization of the initial reactant state (R state) is first performed with a starting configuration chosen from MM/MD-generated structures for the pathways with five water molecules at the dihedral angle of ∼100°. Two sequential paths, C2 and C3 (defined

|

(2)

is then determined. The potential energy profile and the optimized structures of the QM region along this reaction path are shown in Figure 7a and b-d, respectively. The proton is transferred from an oxygen atom of Asp-L210 and reattaches to the same oxygen atom by migrating through a hydrogenbond network around Asp-L210. This proton transfer causes the rearrangement of the hydrogen-bond network of W5-W3 and the orientation of carboxyl group of Asp-L210. The formation of this intermediate state (I state) is almost equienergetic (∆E ) 1.3 kcal/mol) with the R state and yields a moderate energy barrier (6.8 kcal/mol). There exists a very strong HB between Asp-L210 and a water molecule (W3) also in the I state (Figure 7d), giving rise to a very short O-O distance 2.50 Å, even shorter than that of the R state (Figure 7b). This extremely short O-O distance can be attributed to an electronic structure resonance between AspL210-H-H2O and Asp-L210--H3O+, which induces a large low-frequency shift of the O-H vibration.36,41,50-53 Note that the latter electronic structure, Asp-L210--H3O+, is stabilized with the electronically negative environment of the acidic residues such as Asp-M17, Asp-L213, and Glu-H173 and may become a precursor of the formation of a protonated water cluster in the proton channel.6,27,54 Such a resonance structure may also exist in the QAQB- state and explain why FTIR experiments6,29-31 did not detect a distinct peak of the protonated carboxyl group in the proton channel. Near the transition state along this path, a hydronium ion transiently forms HB with Asp-L213. A negative charge transfer then takes place from Asp-L213 to the hydronium ion, which significantly weakens the electrostatic interactions of Asp-L213 with the surrounding positive charged resides, especially ArgL217. These interactions, appearing as the part of the QM/MM interaction, mainly contribute to make the energy barrier at the transition state, as seen in Figure 7a. The proton uptake of QBH- is completed by a subsequent proton transfer starting from the I structure. A minimum-energy path for this proton-transfer process is determined along a reaction coordinate, C3, defined as

C3 )

(| |





| |





RO1 - RH1 - RO2 - RH1





RO3 - RH2

|) (| +





|) (| | | +





|

RO2 - RH2 -





RO4 - RH3 - RO5 - RH3

|)

(3)

The potential energy profile along this proton-transfer path is plotted in Figure 8a. The optimized structures of the QM region are depicted in Figure 8b-d. As seen in Figure 8a, the proton transfer from Asp-L210 to QBH- is exothermic (∆E ) -4.9 kcal/mol), and the activation barrier at the transition state (T2) is 11.9 kcal/mol. This transition state (T2) was determined

Reactions in Reaction Center of Rhodobacter sphaeroides

J. Phys. Chem. B, Vol. 113, No. 26, 2009 8999

Figure 7. Proton transfer through the pathway with five water molecules at 100°. (a) The potential energy profile along this pathway. Blue, red, green, and cyan lines indicate the energies of the overall system, QM, QM/MM interaction, and MM, respectively. The R, T1, and I indicate the reactant, the transition, and the intermediate states, respectively. (b) The QM/MM-optimized structure of the R state. The defined QM region containing the ring part of the QBH-, the side chain of Asp-L210 and Asp-L213, and five water molecules is depicted. The QM/MM-optimized structures of the QM region at T1 and at I are shown in (c) and (d), respectively.

Figure 8. Proton transfer from the intermediate state (I) to the product state (P). (a) The potential energy profile along this path. Blue, red, green, and cyan lines indicate the energy of the overall system, the energy of QM, the QM/MM interaction, and the energy of MM, respectively. I, T2, and P indicate the intermediate, the transition, and the product states, respectively. The QM/MM-optimized structures of the QM region at I, at T2, and at P are shown in (b), (c), and (d), respectively.

by a saddle point search, starting from a potential energy barrier top of the C3 minimum-energy path. The frequency of the imaginary mode at T2 is found to be 565i cm-1. Assuming the reaction coordinate is the direction of O-H stretching, with an “attempt frequency”, and the harmonic approximation of the reactant-state potential, we use a simple transition-state theory to estimate the rate. The attempt frequency characterizes how often the system vibrates and attempts to leave the reactant potential well. The frequency ν of the H-O stretching of a water

molecule in water is 3300-3500 cm-1. Using this wavenumber (∼3300 cm-1) and the value of the speed of light (∼3.0 × 10-2 cm/ps), the frequency (the prefactor) is calculated to be ∼99.0 ps-1. A transition-state theory using this activation energy barrier height and assuming a prefactor of 100 ps-1 predicts a reaction rate to be on a microsecond time scale, which is consistent with experimental observations.49 The potential energy increases as the proton detaches from Asp-L210 with the strong proton affinity. At the transition state

9000

J. Phys. Chem. B, Vol. 113, No. 26, 2009

Kaneko et al.

Figure 9. Relaxation process after the proton transfer. (a) The potential energy profile along the path. Blue, red, green, and cyan lines indicate the energies of the overall system, QM, QM/MM interaction, and MM, respectively. P, T3, and RP indicate the product state, the transition state, and the relaxed product state, respectively. (b) The QM/MM-optimized structure of the P state. The defined QM region containing the ring part of the QBH-, the side chain of Asp-L210 and Asp-L213, and six water molecules is depicted. The QM/MM structures of the QM region at T3 and RP are shown in (c) and (d), respectively.

(T2), a hydronium ion at W2 makes the Eigen structure, rather than the Zundel-structure, with surrounding water molecules in the negatively charged region surrounded by Asp-L210, AspL213, QBH-, and Glu-H173. The hydronium ion attracts a partial electron from Asp-L213, making the latter electronic negativity small. This results in a weakening of the electrostatic interaction of Asp-L213 with surrounding positive charged resides, especially Arg-L217, which is seen in the QM/MM interaction energy of Figure 8a. This destabilization is even enhanced right after T2, that is, at around C3 ) 0.78, when the hydronium ion locates at W1 and makes the direct hydrogen bond to Asp-L213. Here, the partial electron transfer from QBH- to the hydronium ion largely decreases electrostatic stabilization of QBH- by surrounding residues, mainly the Fe2+ ion (data not shown), which is located between QA and QB (this electrostatic interaction appears as the QM/MM interaction in Figure 8a). As the proton further approaches QBH-, the potential energy greatly decreases due to a strong proton affinity of QBH- (the QM energy in Figure 8a). At the product (P) state, Asp-L213 recovers a full electron charge and regains the strong interactions with surrounding positive charged group, particularly with Arg-L217. The QM/MM interaction, therefore, again stabilizes the P state. It should be noted that, although the sequential proton-transfer paths C2 (Figure 7) and C3 (Figure 8) were determined at the ab initio QM/MM level of theory, the starting structure characterized in the distribution of the HB paths shown in Figure 4 was taken from the MD simulation for the protonated AspL210 state with MM force fields, which limits accuracy of the HB paths. In order to find more exact paths, MD simulations using a QM/MM description, instead of MM/MD simulations, shall be performed along the proton transfer, which is beyond the scope of the present calculation.55 Actually, there exists a proton-transfer path (not shown in figures) directly going through W1-W5 without having an intermediate state, I. In this direct path, W3 is twisted in some

angle so that the proton transfer from W3 to Asp-L210, shown in Figure 8c,d is closed, and then, an excess proton moves to W2-W1 directly but yields a much higher total energy barrier than that involved in C2 and C3 paths. Spectroscopic experiments for a mutant, E(L212)Q, have shown that QBH- is accumulated even in the absence of a proton donor at the L212 position,16,38,56 suggesting that there exist alternative proton-transfer pathways other than C1 for the first proton uptake of QB. The transient pathways of C2 and C3 for the second proton uptake of QB may serve as such a path in E(L212)Q. This is examined by constructing the I state for the mutant E(L212)Q and determining the proton-transfer path (the C3 path) from this optimized structure before and after the second electron transfer. The QM/MM-optimized structure in the QA-QB- state has a very weak hydrogen bond between W1 and QB. The proton transfer through this weak hydrogen bond in the QA-QB- state yields the monotonic energy increase of ∼18 kcal/mol. On the other hand, the activation barrier of the proton transfer in QAQB2- is determined to be ∼10.2 kcal/mol, and the reaction is found to be exothermic (∆E ) ∼-9.8 kcal/ mol). With this value of the activation energy, the reaction time for this transfer from the I state is estimated to be on a microsecond order. In the mutation of S(L223)A, the hydrogen bond between W1 and quinone is found to be broken in the MD simulation, and thus, the proton-transfer pathway cannot be formed in the proton channel. The decay time of the semiquinone (QB-) absorption spectrum was found to be several hundred microseconds. This decay time is considered to be the time required for the second electron transfer and the concomitant proton transfer to QB-. In order to calculate such a decay time, we must evaluate the proton-transfer rate by explicitly considering its coupling with the electrontransfer process, but it is beyond the scope of the present work. The P state obtained by QM/MM optimization is further stabilized by conformational relaxation mainly taking place

Reactions in Reaction Center of Rhodobacter sphaeroides

J. Phys. Chem. B, Vol. 113, No. 26, 2009 9001

Figure 10. Proton transfer through the pathway with seven water molecules at 180°. (a) The potential energy profile along the path. Blue, red, green, and cyan lines indicate the energies of the overall system, QM, QM/MM interaction, and MM, respectively. (b-d) The QM/MM-optimized structures along the reaction coordinate. The defined QM region containing the ring part of the QBH-, the side chains of Asp-L210 and Asp-L213, and seven water molecules is depicted.

Figure 11. Schematic diagram of potential energies along the protontransfer path with five water molecules from Asp-L210 to QBH-. R, T1, I, T2, P, T3, and RP indicate the reactant state, the first transition state, the intermediate state, the second transition state, the product state, the third transition state, and the relaxed product state shown in Figures 7-9.

and yields a very small activation (1.5 kcal/mol); the facile conformational relaxation must take place. The proton-transfer pathway with five water molecules is also identified at the Asp-L210 dihedral angle of ∼-75° in Figure 4c. The structural changes associated with this HB pathway at -75° are almost identical to those at 100°, except for the conformational change around the protonation sites of the carboxyl group in Asp-L210. The potential energy profile of the proton transfers at -75° (with the energy barrier height of ∼11.4 kcal/mol, estimated from QM/MM calculation) is similar to those at 100° shown in Figures 7-9. A proton-transfer reaction profile along the longer pathway with seven water molecules at 180°, where the distribution of the paths exhibits its maximum (see Figure 4c), is examined by the QM/MM calculations. Figure 10 depicts the reaction profile and the optimized structures of the initial, transition, and final states. A minimum-energy path is computed along a reaction coordinate, C5, defined as

C5 ) around Asp-L210. A MM/MD trajectory simulation at a very low temperature (10 K) is carried out from the P state configuration in order to escape from a shallow local minimum of the P state. Within 5 ps of the MM/MD simulation, AspL210 undergoes the orientational relaxation, resulting in a new HB formation between Asp-L210 and a water molecule. The minimum-energy path from the P state to the relaxed product (RP) state is computed with the QM/MM method along a reaction coordinate, C4, defined as

|





C4 ) RO1 - RO2

|

(4)

The potential energy profile and the QM region structures in this process are shown in Figure 9. This conformational change from P to RP is found to be exothermic (∆E ) -5.1 kcal/mol)

(|





| |





RO1 - RH2 - RO1 - RH1

|) (| | +





|

RO2 - RH4 -





RO2 - RH3

|)

(5)

One can see in Figure 10 that the reaction process along this pathway requires large energy activation, more than 20 kcal/ mol, indicating its considerably slower reaction rate in comparison with that along the pathway with five water molecules described above. Other sets of reaction coordinates are also examined, but all of the activation energies are found to be larger than 20 kcal/mol. Hence, it can be concluded that, although the formations of the pathways with seven water molecules are more frequent than those with five water molecules, the proton transfer does not proceed through the pathways with seven water molecules. The entire energy profile of the second proton-transfer process is summarized in Figure 11. The rate-determining step is the

9002

J. Phys. Chem. B, Vol. 113, No. 26, 2009

barrier crossing at the second transition state, T2, with the activation energy of 11.9 kcal/mol. In addition, there exists a free-energy cost for the formation of the proton-transferable pathway, which is estimated to be ∼4 kcal/mol based on the relative population of the pathways summed over 100-130°, ∼0.0016, calculated from Figure 4. The overall activation energy of the whole proton-transfer process is therefore ∼16 kcal/mol, which gives a reaction rate on the milliseconds scale, consistent with experimental data.49 4. Conclusions We have investigated the mechanism of the proton transfer in the reaction center of the photosynthetic bacteria. Locations of the excess protons were identified to be at Glu-L212 and Asp-L210 in the QA-QB-state. The mechanisms of the sequential two-proton transfers from these amino acid residues, taking place when an electron is transferred from QA-QB- to QAQB2-, were investigated by using MM/MD calculation and QM/MM optimization. The first uptake mechanism of the proton by QB from Glu-L212 is relatively simple; the proton was found to transfer on the HB directly formed between Glu-L212 and QB. Since the shorter pathway C1 has no potential barrier after the second electron transfer, it is kinetically probable that the shorter pathway is coupled with the second electron transfer in the wild type. This direct proton-transfer path could be altered if GluL212 is highly hydrated. In such a case, there is expected to be a certain energy barrier between these two configurations. Whether the system takes the former or latter configuration must depend on how the proton uptake takes place from bulk water to this moiety upon the first electron transfer. The second proton transfer was found to involve much more complicated steps. The second proton attached to Aps-L210 in the R form interchanges to the I form through the transient HB bridge consisting of three water molecules. It was found that the proton then undergoes the long-range transfer on the HB bridges connecting between Aps-L210 and QB, in the moiety of the proton channel. Formation of such long-range HB bridges occurs very seldom, with a probability of less than 1.0 × 10-3 of the sampled configurations in the trajectory. The barrier energy heights of the proton transfer along these paths, calculated by the QM/MM method, were found to be in good agreement with experimental results, supporting the mechanism proposed in the present study. Further studies are, however, required for understanding the precise mechanism of photochemical reactions in the reaction center. Unless the structure and dynamics of the hydronium ion are well-described, the proton-transfer process along the HB network, consisting of water molecules and amino residues, will not be fully determined. It will be necessary to perform the thermodynamical integral calculation for the free-energy evaluation and the molecular dynamics calculation by using the high-level QM description to treat the hydronium ion properly.57,58 One also needs to explore the proton-transfer mechanism on the multielectronic surfaces for the system with a large QM region, containing QA, QB, several amino acid residues, and water molecules, in order to explicitly examine the coupling between the second electron transfer to QB and the protontransfer dynamics and to determine the total reaction rate of the proton uptake in the photochemical processes of the reaction center.59,60 Furthermore, although the present free-energy calculation has indicated that the protonated R form of Asp-L210 is more favorable as an initial state of the second proton donor, the second proton migrating from bulk water, however, may attach directly to Asp-L210 in the I form rather than that in the

Kaneko et al. R from. In order to identify the initial protonated structure, MD calculations for proton-migration processes from bulk water to Asp-L210 should be performed.61 Acknowledgment. This study was supported by Grants-inAid for Scientific Research (No. 19350009). S.H. was supported by CREST, Japan Science and Technology Agency, and grants from the Ministry of Education, Culture, Sports, Science and Technology of Japan (No.18074004). This study was partly carried out by using the supercomputers at the Research Center for Computational Science in Okazaki, Japan. Y.K. thanks the financial support from the G-COE program of the Nagoya University Chemistry Department. Supporting Information Available: A picture of the movable region in MM/MD simulations. This material is available free of charge via the Internet at http://pubs.acs.org. References and Notes (1) Okamura, M. Y.; Feher, G. Annu. ReV. Biochem. 1992, 61, 861. (2) Hu, X.; Damjanvoic, A.; Ritz, T.; Schulten, K. Proc. Natl. Acad. Sci. U.S.A 1998, 95, 5935. (3) Stowell, M. H. B.; McPhillips, T. M.; Rees, D. C.; Soltis, S. M.; Abresch, E.; Feher, G. Science 1997, 276, 812. (4) Koepke, J.; Krammer, E.; Klingen, A. R.; Sebban, P.; Ullmann, G. M.; Fritzsch, G. J. Mol. Biol. 2007, 371, 396. (5) Xu, Q.; Axelrod, H. L.; Abresch, E. C.; Paddock, M. L.; Okamura, M. Y.; Feher, G. Structure 2004, 12, 703. (6) Nabedryk, E.; Breton, J. Biochim. Biophys. Acta 2008, 1777, 1229. (7) Hasegawa, J.; Nakatsuji, H. Chem. Lett. 2005, 34, 1242. (8) Kleinfeld, D.; Okamura, M. Y.; Feher, G. Biochim. Biophys. Acta 1984, 766, 125. (9) Balabin, I. A.; Onuchic, J. N. Science 2000, 290, 114. (10) Xu., H.; Zhang, R.; Ma, S.; Qu, Z.; Zhang, X.; Zhang, Q. Photosynth. Res. 2002, 74, 11. (11) Keller, S.; Beatty, T.; Paddock, M.; Breton, J.; Leibl, W. Biochemistry 2001, 40, 429. (12) Gerencse´r, L.; Maro´ti, P. Biochemistry 2001, 40, 1850. (13) Giachini, L.; Francia, F.; Mallardi, A.; Palazzo, G.; Carpene`, E.; Boscherini, F.; Venturoli, G. Biophys. J. 2005, 88, 2038. (14) Axelrod, H. L.; Abresch, E. C.; Paddock, M. L.; Okamura, M. Y.; Feher, G. Proc. Natl. Acad. Sci. U.S.A. 2000, 97, 1542. ¨ delroth, P.; Paddock, M. L.; Tehrani, A.; Beatty, J. T.; Feher, G.; (15) A Okamura, M. Y. Biochemistry 2001, 40, 14538. (16) Takahashi, E.; Wraight, C. A. Biochemistry 1992, 31, 855. (17) Paddock, M. L.; McPherson, P. H.; Feher, G.; Okamura, M. Y. Proc. Natl. Acad. Sci. U.S.A. 1990, 87, 6803. (18) Takahashi, E.; Wraight, C. A. Biochim. Biophys. Acta 1990, 1020, 107. (19) Paddock, M. L.; Rongey, S. H.; McPherson, P. H.; Juth, A.; Feher, G.; Okamura, M. Y. Biochemistry 1994, 33, 734. (20) Paddock, M. L.; Feher, G.; Okamura, M. Y. Biochemistry 1995, 34, 15742. ¨ delroth, P.; Chang, C.; Abresch, E. C.; Feher, (21) Paddock, M. L.; A G.; Okamura, M. Y. Biochemistry 2001, 40, 6893. (22) Takahashi, E.; Wraight, C. A. Proc. Natl. Acad. Sci. U.S.A. 1996, 93, 2640. (23) Paddock, M. L.; Feher, G.; Okamura, M. Y. FEBS Lett. 2003, 555, 45. (24) Paddock, M. L.; Flores, M.; Isaacson, R.; Chang, C.; Abresch, E. C.; Okamura, M. Y. Biochemistry 2007, 46 (28), 8234. (25) Alexov, E. G.; Gunner, M. R. Biochemistry 1999, 38, 8253. (26) Ishikita, H.; Morra, G.; Knapp, E. W. Biochemistry 2003, 42, 3882. (27) Hermes, S.; Stochnik, J. M.; Onidas, D.; Remy, A.; Hofmann, E.; Gerwert, K. Biochemistry 2006, 45, 13741. (28) Remy, A.; Gerwert, K. Nat. Struct. Biol. 2003, 10, 637. (29) Nabedryk, E.; Breton, J.; Hienerwadel, R.; Fogel, C.; Ma¨ntele, W.; Paddock, M. L.; Okamura, M. Y. Biochemistry 1995, 34, 14722. (30) Nabedryk, E.; Breton, J.; Okamura, M. Y.; Paddock, M. L. Biochemistry 1998, 37, 14457. (31) Nabedryk, E.; Breton, J.; Okamura, M. Y.; Paddock, M. L. Biochemistry 2001, 40, 13826. (32) Cornell, W. D.; Cieplak, P.; Bayly, C. I.; Gould, I. R.; Merz, K. M.; Ferguson, D. M.; Spellmeyer, D. C.; Fox, T.; Caldwell, J. W.; Kollman, P. A. J. Am. Chem. Soc. 1995, 117, 5179. (33) Berendsen, H. J. C.; Grigera, J. R.; Straatsma, T. P. J. Phys. Chem. 1987, 91, 6269.

Reactions in Reaction Center of Rhodobacter sphaeroides (34) Berendsen, H. J. C.; Postma, J. P. M.; van Gunsteren, W. F.; DiNola, A.; Haak, J. R. J. Chem. Phys. 1984, 81, 3684. (35) Ryckaert, J. P.; Ciccotti, G.; Berendsen, H. J. C. J. Comput. Phys. 1977, 23, 327. (36) Hayashi, S.; Ohmine, I. J. Phys. Chem. B 2000, 104 (45), 10678. (37) Schmidt, M. W.; Baldridge, K. K.; Boatz, J. A.; Elbert, S. T.; Gordon, M. S.; Jensen, J. H.; Koseki, S.; Matsunaga, N.; Nguyen, K. A.; Su, S. J.; Windus, T. L.; Dupuis, M.; Montgomery, J. A. J. Comput. Chem. 1993, 14, 1347. (38) Hienerwadel, R.; Grzybek, S.; Fogel, C.; Kreutz, W.; Okamura, M. Y.; Paddock, M. L.; Breton, J.; Nabedryk, E.; Ma¨ntele, W. Biochemistry 1995, 34, 2832. (39) Mezzetti, A.; Nabedryk, E.; Breton, J.; Okamura, M. Y.; Paddock, M. L.; Giacometti, G.; Leibl, W. Biochim. Biophys. Acta 2002, 1553, 320. (40) Brzezinski, P.; Paddock, M. L.; Okamura, M. Y.; Feher, G. Biochi. Biophys. Acta 1997, 1321, 149. (41) Hayashi, S.; Tajkhorshid, E.; Schulten, K. Biophys. J. 2002, 83 (3), 1281. (42) Edman, K.; Nollert, P.; Royant, A.; Belrhali, H.; Pebay-Peyroula, E.; Hajdu, J.; Neutze, R.; Landau, E. M. Nature 1999, 401, 822. (43) Kamiya, M.; Saito, S.; Ohmine, I. J. Phys. Chem. B 2007, 111, 2948. (44) Humphrey, W.; Dake, A.; Schulten, K. J. Mol. Graphics 1996, 14, 33. (45) Verme´glio, A. Biochim. Biophys. Acta 1977, 459 (3), 516.

J. Phys. Chem. B, Vol. 113, No. 26, 2009 9003 (46) Wraight, C. A. Biochim. Biophys. Acta 1977, 459, 525. (47) Graige, M. S.; Feher, G.; Okamura, M. Y. Proc. Natl. Acad. Sci. U.S.A. 1998, 95, 11679. (48) Graige, M. S.; Paddock, M. L.; Bruce, J. M.; Feher, G.; Okamura, M. L. J. Am. Chem. Soc. 1996, 118, 9005. (49) Mezzetti, A.; Leibl, W. Eur. Biophys. J. 2005, 34, 921. (50) Hayashi, S.; Tajkhorshid, E.; Kandori, H.; Schulten, K. J. Am. Chem. Soc. 2004, 126 (34), 10516. (51) Shibata, M.; Kandori, H. Biochemistry 2005, 44 (20), 7406. (52) Shibata, M.; Tanimoto, T.; Kandori, H. J. Am. Chem. Soc. 2003, 125 (44), 13312. (53) Garczarek, F.; Gerwert, K. Nature 2006, 439, 109. (54) Breton, J.; Nabedryk, E. Photosynth. Res. 1998, 55, 301. (55) Riccardi, D.; Kon¨ig, P.; Guo, H.; Cui, Q. Biochemistry 2008, 47, 2369. (56) McPherson, P. H.; Scho¨nfeld, M.; Paddock, M. L.; Okamura, M. Y.; Feher, G. Biochemistry 1994, 33, 1181. (57) Li, G.; Cui, Q. J. Phys. Chem. B 2003, 107 (51), 14521. (58) Ghosh, N.; Cui, Q. J. Phys. Chem. B 2008, 112 (28), 8387. (59) Hammes-Schiffer, S. Acc. Chem. Res. 2001, 34 (4), 273. (60) Borgis, D.; Hynes, J. T. J. Phys. Chem. 1996, 100 (4), 1118. (61) Phataka, P.; Ghoshb, N.; Yub, H.; Cui, Q.; Elstner, M. Proc. Natl. Acad. Sci. U.S.A. 2008, 105, 19672.

JP9008898