M. Macernis a, L. Mourokh b, L. Vinciunas a, and L. Valkunas a,c

aInstitute of Chemical Physics, Faculty of Physics, Vilnius University, Saulėtekio 9-III, 10222 Vilnius, Lithuania

bDepartment of Physics, Queens College of the City University of New York, Flushing, 11367 New York, USA

cDepartment of Molecular Compound Physics, Center for Physical Sciences and Technology, Saulėtekio 3, 10257 Vilnius, Lithuania

Received 22 June 2022; accepted 3 July 2022

Photoinduced proton pumping in bateriorhodopsin (bR) was extensively studied using multiple experimental methods and utilizing various theoretical modelling approaches. These studies usually refer to the well-resolved structural data of bacteriorhodopsin. However, despite the obtained results, the origin of the proton pumping force initiated by the electronic excitation of retinal remains questionable. Using quantum chemical calculations, we have revealed that the retinal molecule after its excitation is fixed in the ground state of 13-cis,15-syn configuration, as a result of interaction with specific protein residuals. Reaching this fixed configuration, the proton is first transferred to the aspartic acid No. 85 (Asp-85) residue from the water molecule, which is subsequently restored by the proton initially located in the Schiff base. We discuss the challenges and approaches to modelling the proton transfer in bR and demonstrate that the process, which starts from the electronic excitation of the retinal molecule, is mainly due to the detailed arrangement of the protein environment.

Keywords: photoinduced proton transfer, membrane proteins, bacteriorhodopsin, retinal, isomerization

PACS: 87.14.ep, 77.65.-j, 82.80.Yc, 83.10.Rs, 03.65.-w

1. Introduction

Bacteriorhodopsin (bR) is the  light-transducing protein in the purple membrane of Halobacterium salinarum [1, 2]. This hydrophobic protein contains seven transmembrane helices (A through G, see Fig. 1), which span the membrane at small angles to the perpendicular and surround the centrally located retinal linked to G-helix by a protonated Schiff base (SB) with Lys-216 [36]. As a  result of the  light absorption, the  retinal yields the  transition into the  13-cis,15-syn form, which triggers a cyclic reaction that comprises the series of spectral intermediates usually named J, K, L, M, N and O, reflecting protein stages of the  proton pumping across the membrane [713]. Proton transport starts from the middle of the protein structure, with the  proton transfer from the  SB to aspartic acid 85 residue labelled as Asp-85 or D85 [14]. Subsequently, the proton is transferred to the extracellular side of the membrane and the SB is reprotonated from the cytoplasmic side completing the cycle. In addition to the  ability of the  proton pumping of bR, it was also experimentally demonstrated that a negative chlorine (Cl) or positive proton (H+) exchanger could be converted into a passive channel for Cl by substitution of only two key residues [15, 16]. Site-directed mutagenesis combined with various spectroscopic methods, crystallographic data, in conjunction with quantum chemical and molecular dynamics simulations, revealed many details of the structural changes of the protein backbone in the course of the bR photocycle [8]. However, the main questions concerning the exact pathway of energy conversion from the  electronic excitation of retinal into the proton pumping force in this photoactive protein, as well as the unidirectionality of proton transfer processes, are still unanswered.


Fig. 1. Schematic presentation of the structure of bacteriorhodopsin: it contains seven transmembrane helices labelled from A through G and retinal linked to G-helix via a protonated Schiff base.

In the present paper, we discuss the challenges for the modelling of bR related to its intermediate size, as it is too large for the fully quantum-chemistry-based description, but nevertheless, the system of helices and residues bound to them cannot be described statistically, in terms of the uniform dielectric environment. Moreover, the energy exchange between various parts of the whole system, such as the  excited electronic states, transferred protons, twisted retinal and bent helices, is not unidirectional, as none of them serves as the energy sink. Correspondingly, the  calculations can be fulfilled by keeping in mind the energetic states of the entire system, even if it is impossible to take into account all the  specific interactions. As an example, we address the  proton transfer event associated with the deprotonation of the SB and the protonation of the Asp-85 residue. We show that the irreversibility of the proton transfer can be achieved in terms of quantum chemical calculations if the 13-cis,15-syn configuration of the retinal defines the initial conditions for the proton pumping. This isomerization does not correspond to the  minimum of energy for the retinal separated from the rest of the complex. However, we argue that it might be defined as the state corresponding to the proper energy redistribution in the full system at this configuration when all interactions with the residues at the bent helices are taken into account. We would also like to stress one of the active hypothesis that according to our analysis, proton transport starts from the middle of the protein, with one of the water molecules donating its proton to Asp-85 with a subsequent reprotonation from the SB.

The rest of the paper is structured as follows. In Section 2, we address the challenges of the modelling of the proton transfer. A review of the molecular dynamics and quantum chemical calculations of the photoisomerization is also presented in Section 2. The methods used for quantum chemical calculations are briefly described in Section 3. The  results obtained by modelling the  retinal isomerization and its link to the initial proton transfer steps are presented in Section 4. Discussion of the obtained results is presented in Section 5, and Section 6 contains the conclusions of our work.

2. Challenges in modelling of bacteriorhodopsin

2.1. Charge transfer

Modelling of charge transport is usually based on the assumption that the system can be effectively separated into the dynamical part and the environment. This assumption is at the heart of the Marcus theory of electron transfer with the  transfer rate given by

k ET = 2π | H 12 | 2 1 4πλ k B T exp( ( λ+Δ G 0 ) 4λ k B T 2 ),(1)

where |H12| is the coupling between the initial and final states, λ is the total reorganization energy of the  environment, kB is the  Boltzmann constant, T is the absolute temperature, and ΔG0 is the total Gibbs free energy change for the transfer reaction. It was shown that such an approach can be successfully used to describe proton transport as well [17, 18]. However, there are two caveats for a specific bR case, both of them related to the smallness of the system and the corresponding impossibility of the separation of the dynamical subsystem and its environment. First, it is very difficult to define the reorganization energy, with all the bR residues reacting to the  proton transfer events in a  rather dynamical, not statistical manner. Second, in the  calculation of the  differences of free energies between the initial and final states, the direct electrostatic interaction between the transferring proton and charges on the residues provides an essential impact on the  result. In previous works based on the energy conceptions [19, 20], these both aspects are not taken into account properly. Free energies corresponding to proton sites (based on the related pKa factors) are usually calculated with the  electrostatic interactions considered. The rest of the system is included by calculating the  reorganization energy and non-vacuum dielectric permittivity, the  electrostatic interaction of the protons with electrons being completely ignored, in particular. Water molecules, which are shown to be crucial for the proton transfer in bR [2123], are considered electrically neutral, non-dissolved and dipole-less. While this approach can be appropriate for large protein systems [24], its applicability for relatively small systems, such as bR, is questionable.

A proper description of the  full bR dynamics could be based on the inclusion of the influence of all helixes and residues into the  Gibbs free energy and the  consideration of various vibrational degrees of freedom as the environment. In this case, the reorganization energy would be small, and the functional form of Eq. (1) becomes close to a delta-function. With small dissipation, the  transfer between the  states would be almost resonant, but slightly downhill with respect to the energy of the full system. As mentioned above, the bR complex is too large for such quantitative treatment, but this approach might be used for qualitative discussion. In particular, if the  exact dynamics of some fragment of the  full complex is being calculated, the  initial and final energies of this fragment alone do not have to be minimal for given configurations, as the  requirement of the energy minimum for the relaxed system applies only to the whole complex. Moreover, when two final states are available, the system can have a higher probability to be transferred to the more energetic of the  two, because it is impossible to unload a large amount of energy to the vibrational environment in order to proceed to the minimal energy state, i.e. if the  system is in the  inverted Marcus regime. In some sense, our approach is comparable to the  ‘reaction coordinate models’ [25, 26] where strongly-coupled environmental degrees of freedom (often high frequency modes) were specifically included into the  dynamical part.

2.2. Molecular dynamics and quantum chemical calculations

Mimicking photoinduced processes in bR is not a trivial task since the resulting outcome is sensitive to the protein environment, that has to be taken into account properly. For this purpose, molecular dynamics (MD) simulations have been intensively used [4, 2737]. Indeed, the rotation of SB as a  result of the  retinal isomerization (see Fig. 2) can be obtained by modelling the  retinal with its protein environment [34, 35]. However, in order to achieve it, multiple approximations should be used due to the limits of computational resources including algorithms themselves. One very crucial approximation suggests that the retinal structure can be subdivided into quantum mechanics (QM) and molecular mechanics (MM) segments, and the division between them is set across covalent bonds. For QM/ MM simulations various approximations can be used for the  description of MM and QM counterparts. The  definition of the  QM–MM boundary is the major problem by constructing covalent bonding and adaptative schemes [4, 27, 33, 37]. Hayashi and co-authors [34, 35] exploited only a part of the retinal structure with the attached SB using the QM approach, while the rest of the molecule was parameterized. QM simulations of full retinal do not demonstrate isomerization and/or rotations of some parts of the molecule as follows from experimental data [36, 37]. Moreover, the SB bound to the retinal does not rotate in the vacuum models but instead only the rotation of the β-ring is obtained in the excited state of retinal [38, 39]. Starting from the excited state, the structure optimization procedure restores the  initially optimized ground state structure, which implies that the  minima of the  potential energy surfaces in the ground- and excited-states are not separated by any significant energy barrier. Furthermore, by analyzing the isomerization process of retinal up to 90 deg, it was concluded that the DFT methodology is not applicable due to the fact that the system is reaching the conical intersection region in the course of isomerization [40].

Molecular dynamics simulations with molecular mechanics force fields are performed in order to relate deformations of the bR helices to the proton transfer mechanism by assuming that the  L-intermediate follows from the  retinal alltrans to 13-cis photoisomerization step corresponding to the K-intermediate [41, 42]. If the water molecules (especially W402) are not taken into account, the retinal pumping of the proton is expected from the 13-cis,15-anti structure and not from the 13-cis,15-syn, in opposite to the crystallographic data [42]. Other studies based on three transfer mechanisms in terms of QM/MM calculations [43] also indicate the presence of the 13-cis,15-anti configuration. More accurate approach for the simulation of the protonated retinal active centre is based on the force field methodology. The water molecular clusters with additional possible water molecules in accordance with several crystallographic data were taken into consideration [44]. As follows from this type of studies, the proton transfer from retinal interacting with the water molecule (W402) is more complicated since the presence of the energy barrier appears for the proton pathway [45], while the NMR studies report opposite evidences [46].


Fig. 2. Schematic representation of the photoinduced initial processes: initial state of the retinal (RET) coupled with the protonated Schiff base (SBH+); photoinduced isomerization of the retinal and subsequent proton release [5, 40]. The Schiff base (SB) rotation angle presents the dihedral angle described by atoms 12-13-14-15.

The protonated retinal after its excitation is known as the J state of bR corresponding to a complication of the excitation dynamics of the retinal [45-47]. In the  vicinity of the local minimum of the first singlet excite state (S1), the retinal coupled to the SB remains almost planar, but the single and double bonds turn to become inverted according to the  MD simulations and ab  initio modelling [29, 30, 32, 37, 48]. As follows from these data, the  existence of a  small barrier between trans and cis structures in the excited states could partially explain the fluorescence effects related to the structural evolution of retinal bound to the  protonated SB in bR [29]. Our studies of the retinal-type conjugated systems support this conclusion: the  excited state S1 demonstrates an inversion between single and double bonds, that should be caused by the  presence of the  second singlet excited state (S2) properties (also known as 11Bu+ type state, corresponding to its symmetry according to polyene models), appropriate for polyene and carotenoid molecules [49]. It follows from the equation-of-motion coupled cluster singles and doubles model (EOM-CCSD) approach [50] that the  lowest excited state of the  retinal bound to the protonated SB is of a similar origin to that of the polyene type 11Bu+ state. In the case when retinal is bound to the non-protonated SB, the  optically forbidden state of symmetry 11Ag turns to be the lowest excited state due to the electron correlations. The  EOM-CCSD is computationally expensive and the excited state energies are largely overestimated. The  lowest S1 state of the  retinal bound with the  deprotonated SB should be of the 11Ag symmetry while the lowest S1 state of the retinal bound to the protonated SB should be of the  11Bu+ symmetry. Thus, the  retinal energy surfaces should be more complicated by considering the retinal isomerization together with the  possible flexibility of helices. Thus, energy transfer channels, which might fix the retinal in the 13-cis,15-syn and 13-cis,15-anti configuration (Fig. 2), and finally pass the absorbed light energy from the  retinal into the  helices, are expected.

Deprotonation of the  SB is usually attributed to the first proton transfer step in bR [10]. The protonation of Asp-85 upon the deprotonation of SB together with the photoisomerization of retinal from the initial all-trans configuration to the 13-cis configuration is followed by the release of a proton via a network of protein residues and bound water molecules [5156] to the extracellular membrane surface [57, 58]. As follows from studies based on theoretical calculations, various possible pathways for retinal deprotonation in the ground state, which at moderate approximations agrees with the  experimental data [48, 59], are expected.

3. Methods

To reveal the connection between the retinal photoisomerization and the  initial steps of proton transfer, the studies of possible energy redistribution channels between the electronic excitation of the  retinal and the  protein surrounding are carried out. For these purposes, stationary density functional theory (DFT) calculation, two-layer own N-layered integrated molecular orbital and molecular mechanics (ONIOM), and quantum molecular dynamics (QMD) are used. Recently it was demonstrated that by considering simple polyene type molecules the excited states of retinal have a multireference character, which cannot be taken into consideration within the  DFT approach [45, 60, 61]. However, the proton transfer is taking place from the ground state of the retinal as a result of the excitation relaxation, thus causing several bR stages [1, 9, 18, 53, 62, 63]. Therefore, the DFT approaches should be applicable to provide quite accurate descriptions.

It is known that the bR in the L state reflects changes of the active centre environment in comparison to the initial state [62]. Thus, to simplify the description of the initial steps of proton pumping, these two states will be mainly considered. The initial structures are taken from the protein data bank (PDB). The  PDB codes for the  structures are 1C3W [5] and 2NTW [64], respectively. The  structure 1C3W provides crystallographic data for bR in the  initial state with the  all-trans configuration of the retinal and the 2NTW provides the crystallographic information for the bR in the  intermediate state with the  13-cis,15-syn configuration of the retinal (Fig. 2).

For QM calculations, DFT modelling is performed using the Gaussian software package [50] with HSEH1PBE, CAM-B3LYP and B3LYP functionals, and cc-pVDZ and 6-31G* basis sets. Similar approaches were used previously [27, 46, 49]. It should be noted that the HSEH1PBE and CAMB3LYP functionals involve long-range corrections. According to our studies, the  HSEH1PBE functional provides the same results as the B3LYP functional. Here we report the  results of computations with the  HSEH1PBE functional with the  cc-pVDZ basis set. The  CAM-B3LYP functional, which results in the β-ring rotation when the  all-trans retinal is optimized, is also applied for the  optimization of the  retinal structure in the  excited state. The  QM/MM or QM/QM in terms of ONIOM approach are carried out by adding relevant residuals of the protein and water molecules into the  calculation scheme using molecular mechanics with the  universal force field (UFF) or semiempirically using PM6 Hamiltonian methods. In all cases, the  lower layer was always frozen while the higher layer was free of changes in ONIOM methodology. Quantum molecular dynamics (QMD) is performed with Gaussian 9 version at the  B3LYP/cc-pVDZ level by using the atom-centred density matrix propagation (ADMP) molecular dynamics model at different temperatures: 1, 70 and 300  K. ADMP was performed with 63 atoms taken into account, time steps were 0.1  fs, fictitious electronic mass was 0.1  amu, velocity scaling was performed all the  way through the  simulation, the  dynamics were performed with converged SCF results at each point, and the temperature was set to determine the fixed nuclear kinetic energy according to the Gaussian overlay documentations [50].

4. Results

4.1. Modelling of retinal isomerization

The proton transfer correlates with the photoinduced transition of the  retinal from its initial trans-configuration to the  13-cis,15-syn isomer on µs timescale [65]. Thus, first, we consider the photoisomerization process in terms of QM/MM by taking into account the  structural details of the  part of the  protein from the  vicinity of the photoactive centre (retinal) as follows from crystallographic data [5, 64]. The protonated SB, which is a  proton donor, becomes connected to the  proton acceptor Asp-85 via the  water molecule (labelled W402) as a  result of the  retinal isomerization. Additional water molecules, labelled W401 and W406, as well as the protein two tryptophan residues labelled Trp-86, Trp-182 and tyrosine residue (labelled Tyr-185) are also present in the vicinity of the active centre and might play a decisive role in the photoinduced proton transfer. In order to repeat calculations with a higher basis set the  threonine residue (labelled Thr-90) may also be required. The all-trans configuration of the  retinal together with the  protein residues and bound water molecules mentioned above is taken from the crystallographic data.

As follows from the calculations, the all-trans retinal (corresponding to the  ground state) contains only one global minimum along the rotational pathway, while the optimization in the excited state results in the presence of two possible minima (Fig. 3(a)). One of these minima is at about –92 deg that correlates well with the structure defined from the  crystallographic data 2NTW [64]. Moreover, the ground state in this configuration can also be found when the protein residues Trp-86, Trp-182, and Tyr-185 are taken into account, thus, suggesting a conical intersection between these two states (Fig. 3(b)). Similar results are obtained both in QM/MM or QM/QM (ONIOM) studies by adding Trp-86, Trp-182, and Tyr-185 residues into the calculation scheme and using molecular mechanics UFF and in the semiempirical PM6 calculations. Thus, the minima corresponding to the excited states of 13-cis,15-syn isomer of the retinal appear within the supermolecular (several molecules in the same QM calculations) approach containing the retinal Trp-86, Trp-182, and Tyr-185 residues in agreement with the earlier suggestions for the ground state stability. However, this structure does not follow from the calculations of the ground state, so it can be obtained for the excited state minima only (Fig. 3). Accordingly, to obtain the minimum for the 13-cis,15-syn isomer of retinal in the ground state, additional interactions should be taken into account.


Fig. 3. Schiff base rotation energy surface of the retinal bound to the protonated Schiff base along the dihedral angel indicated in Fig.  2: 0 deg corresponds to the all-trans retinal configuration in vacuum (a); –92.4 deg represents the 13-cis,15-syn retinal isomer together with the Trp-86, Trp-182 and Tyr-185 residues in vacuum (b). The dots represent the dihedral scan increment.

Thus, to examine the stability of the configuration of active centre corresponding to the retinal in the 13-cis,15-syn configuration, we initially fix this configuration and use the QMD simulations of the structure at different temperatures (Fig. 4).

The initial structures used for the calculations are chosen to be optimized either in the ground state or in the first excited state with the SB rotated by 90 deg. The results presented in Fig. 4 demonstrate that the ground state structure (E0) remains almost the same, while in the excited state (E1), this structure reaches other values of the rotational angle (retinal of cis isomer) on the time scale less than 0.5  ps at low temperatures. At higher temperatures (70 and 300  K), the  isomerization process starts at around 30  fs and experiences transition into another structural arrangement in about 100 fs and, thus, the initial state evidently cannot survive for several µs, as expected from the  experimental observations. Accordingly, we conclude that QMD simulations do not support the  evidence of the  appearance of 13-cis,15-anti retinal isomer without taking into account the essential interaction with the protein environment that has been ignored so far.

4.2. Initial step of proton transfer

To model the initial stage of the proton transfer, let us keep the key atoms taken from the crystallography data (with the  retinal in the  13-cis, 15-syn configuration and water molecules as it is in 2NTW [64]) fixed (see Fig. 5). The modelling data obtained by using QM calculations demonstrate the  proton displacement from W402 to the Asp85 group (Fig. 6(a)) that opens the pathway of the proton transfer from the SB to W402 (Fig. 6(b)).


Fig. 4. Possible evolution of the retinal demonstrating the isomerization as follows from the ADMP type quantum molecular dynamics in the ground state (E0) and in the first excited state (E1) at 1, 70 and 300 K temperatures, respectively. Initially, the optimized 13-cis ,15-syn configuration is assumed (Fig. 3(b)). Here α denotes the rotation angle (see Fig. 2). The inset shows the simulation results for the structure in the E1 state at 1 K up to 1 ps.


Fig. 5. Active centre, protein residuals and bound water molecules which are taken into account for the calculations. The retinal bound to the protonated Schiff base is rotated by –92.4 deg of the dihedral angle (see Fig. 2), in agreement with the structural data of the state of bR when the proton transfer occurs.

Since the proton transfer is initiated by the photoexcitation of the retinal bound to SB, it subsequently rotates towards the  active centre where the water molecule mediating the Asp-85 residue is positioned, as shown in Fig. 5.  Thus, the photoexcitation of the retinal molecule causes the transition from the trans retinal (0 deg) to 13-cis configuration (–92.4 deg) in the excited S1 state.

In this case, the  proton transfer starts from the  water molecule (W402) to Asp-85 residue with a subsequent proton transfer from the protonated SB to the  formed temporary hydroxide (OH). This explains the  pathway demonstrating the  role of water molecules as follows from the energy surface scans shown in Fig. 6(a). It should be mentioned that the presence of water molecules W401 and W406 is crucial for the proton transfer from the W402 to the Asp-85 group (Fig. 6(a)). These molecules make an essential influence on the energy surfaces by opening protonation pathways in bR.

5. Discussion

5.1. Energy redistribution during the photocycle

Let us discuss now how the  energy supplied by the  photon is redistributed over various degrees of freedom in the  transitions between the  intermediates shown in Fig.  7, keeping in mind that these transitions are almost resonant with a small dissipation. Immediately after the photon absorption, the  energy is entirely in the  electronic degrees of freedom of the excited state of the retinal. The excited electronic state cannot relax back to the all-trans ground state because the internal degrees of freedom are not able to accept all the energy of the electronic excitation. Instead, the conformation to the  13-cis form (J to K transition) allows the excitation energy to be transferred to the  strain of the  retinal. It should be noted that the part of the energy is still contained in the electronic subsystem, as the  ground state energy of the 13-cis isomer is larger than that of all-trans. In the consequent transition between K and L intermediates, the retinal is more relaxed and the energy is redistributed to the  protein resulting in bending of the C-helix and thus in modification of the energies of the Coulomb interaction between residues [63, 66].


Fig. 6. The proton pathway from the protonated SB (in the 13-cis, 15-syn configuration of retinal) to Asp-85 via mediating water molecule (W402): (a) energy surfaces between the water molecule (W402) and Asp-85 at the presence of the 13-cis, 15-syn retinal and other water molecules (W401, W406); (b) energy surfaces between the SB and water molecule (W402) due to the presence of Asp-85 group. Here RN-H defines the distance between the nitrogen atom in SB and the hydrogen atom, as shown in Fig. 5. The computations were performed using the HSEH1PBE functional.


Fig. 7. The  photocycle of bR showing the  main proton transfer events and the main mechanical transformations of the protein and its environment.

These modified interactions compensate for the  huge acid dissociation constant (pKa) differences between the SB and Asp-85 (more than 11 and 2.5, respectively, for the bR states) and allow the  initial proton transfer between these residues (L to M1 transition). To make the transfer possible, the pKa factors of both protonable sites are modified bringing them to the pKa of about 6.5. Accordingly, the energy delivered to the proton is much less than half of the electron-volt as follows from the initial pKa difference. The protonation of Asp-85 modifies the electrostatic interactions in the proton release group, so the  corresponding proton energy becomes larger than the chemical potential of the extracellular space and the proton can be released (M1 to M2 transition). The deprotonation of SB leads to conformational changes (relaxation of the C-helix and bending of B- and F-helices), which, in turn, allow both the increase of the pKa of SB to 8.3 and the decrease of the pKa of Asp-96 to 7 [8]. Moreover, it leads to the penetration of water molecules in the space between the Asp-96 residue and SB. As a result, the reprotonation of the SB from Asp-96 via the water bridge becomes energetically favourable (M2 to N1 transition). Consequently, the proton energy level of Asp-96 again becomes smaller than the chemical potential of the cytoplasmic side and the Asp-96 residue is reprotonated (N1 to N2 transition). Accordingly, in the  series of events the  bR complex releases the  high-energy proton to the extracellular side and takes the low-energy proton from the  cytoplasmic side, i.e. the  proton is effectively pumped. The energy gain can be estimated from the difference between the pKa of Asp-96 residue at the initial state (about 11) and the potential of the hydrogen (pH) factor of extracellular space, which is about 0.3 eV for pH = 5. The proton translocations from and to Asp-96 give rise to a new sequence of conformational changes, such as the relaxation of B- and F-helices, and the modifications of electrostatic interactions, so the retinal reisomerization to the  all-trans configuration is now possible (N2 to O transition). The reprotonation of SB leads to the decrease of the pKa of Asp-85 and its proton is slowly transferred downhill to the proton release group (O to bR transition) completing the photocycle.

5.2. General description of energy redistribution processes

Let us discuss now how the obtained results of our calculations of the initial step of proton transfer are related to a  general description of the  energy redistribution processes. Proton transport correlates with the protein conformational changes at various bR stages during the  photocycle, especially with the  bending of C- and F-helices (Fig.  7). In particular, the proton transfer to the extracellular side is related to the C-helix, while the proton uptake from the cytoplasmic side is related to the F-helix.

As follows from our quantum chemistry studies, the rotation energy surface cannot explain the initial energy transfer without the account of interaction with the protein helices. Indeed, the protonated 13-cis,15-syn retinal configuration cannot occur without the support of the protein because it is absent in the ground state structure (Fig. 3(a)), while the relevant minima might be present in the first excited state if three structural parts from the bR are taken into account (Fig. 3(b)). Thus, the energy should be evidently translocated to the C-helix for the rearranged positions of three water molecules (Fig. 6(a)), which appear as a result of the isomerization in the vicinity of the protonated Schiff base. This complex dynamics of molecules and protein residuals in the vicinity of the excited retinal could result in irreversible proton transfer by (i) temporary shifting the  OH- group from W402 and subsequently protonating Asp-85 (Fig.  6(b)) and (ii) translocating the proton from the SB group to W402 (Fig. 6(a)). Thus, the initial energy should be transferred to the protein in the initial step.

As a  result of such configuration changes, the energy of the electronic states of retinal should be modified as schematically shown in Fig. 8. According to our calculations, the proton should be transferred to Asp-85 (with the help of the water molecule) from the  13-cis,15-syn configuration, which is not of the lowest energy state of the retinal taken alone. To have the retinal fixed in this configuration, we have to postulate that the energy of this state is decreased as a result of the electrostatic interaction with the bent helices. The ground state of the  retinal at the  13-cis,15-syn conformation with the account of this interaction should be almost at the same energy level as it is for cis/trans minima. Thus, we should have two similar centres, which differ by the  water translocation but with similar structural parts causing different irreversible proton transfer events (Fig. 8).


Fig. 8. Energy model for the retinal electronic excited states and possible excitation energy uptake by the protein of bR. The red lines mark the energy surfaces of the retinal, the green lines correspond to the energy changes of the protein.

The energy translocation from the retinal to protein helices in the primary step allows us to qualitatively determine the energy swap in the protein, as shown in Fig. 8. Proton transport is initiated from the middle of the structure by step-by-step pushing and uptaking protons over all the bR proteins resulting in pKa changes.

6. Conclusions

In conclusion, we have discussed the challenges and approaches to modelling the proton transfer in bR and suggested that the retinal molecule responsible for the light absorption cannot be completely separated from the protein environment. In the spirit of the  Marcus transport theory, all the  energies of the initial and final states of the transfer event should be taken from the whole system, including the  energies of the  bent helices, stretched bonds and modified electrostatic couplings to the charged and electric-dipole-possessed residues and water molecules. If the  part of the  whole system is involved in specific exact or high-level numerical calculations, it does not have to be at its ground state (calculated for this isolated part), as one of the excited states can correspond to the energy minimum for the whole system.

For the illustration, we have discussed molecular dynamics and performed quantum chemistry calculations of the initial proton transfer event in bR, i.e. the deprotonation of SB and the protonation of Asp85 residue, with the retinal isomerization and the energy of the  retinal electron state taken into account. We have shown that the proton transport becomes irreversible if the retinal is initially in the  13-cis,15-syn configuration, which is not the  ground state. However, this state might correspond to the energy minimum of the whole system including the energy of the modified electrostatic interaction with the  residues of the  bent C-helix. We have revealed the actual mechanism of the proton transfer, as initially the Asp-85 residue is protonated by the water molecule which is subsequently restored by the SB proton. In some sense, it is comparable to the  proton transfer through the whole complex, as it is not the situation when the single proton is transferred from the cytoplasmic side to the  extracellular space, but instead, proton transport is initiated from the  middle of the structure.

We believe that the outlined approach will offer additional prospects for future research because the restrictions on the initial and final states can be lifted. Correspondingly, molecular dynamics and quantum chemistry calculations can be performed in a much wider range of parameters and allow detailed descriptions of individual transfer events. We also believe that eventually all these individual events will be combined in the unified self-consistent picture of the bR machinery.


This work was supported by the  Research Council of Lithuania (LMT Grant No.  P-MIP-20-219). Computations were performed by using the High-Performance Computing Centre ‘HPC Saulėtekis’ at the Faculty of Physics of Vilnius University. L. M. acknowledges the support by the Air Force Office of Scientific Research, Award No. FA9550-16-1-0279.


[1] D.  Oesterhelt and W.  Stoeckenius, Rhodopsinlike protein from the  purple membrane of Halobacterium halobium, Nat. New. Biol. 233, 149–152 (1971).

[2] R.F. Peck, S. DasSarma, and M.P. Krebs, Homologous gene knockout in the  archaeon Halobacterium salinarum with ura3 as a  counterselectable marker, Mol. Microbiol. 35, 667–676 (2000).

[3] N.  Grigorieff, T.A.  Ceska, K.H.  Downing, J.M.  Baldwin, and R.  Henderson, Electron-crystallographic refinement of the structure of bacteriorhodopsin, J. Mol. Biol. 259, 393–421 (1996).

[4] U.N. Morzan, D.J. Alonso de Armino, N.O. Foglia, F. Ramirez, M.C. Gonzalez Lebrero, D.A. Scherlis, and D.A. Estrin, Spectroscopy in complex environments from QM-MM simulations, Chem. Rev. 118, 4071–4113 (2018).

[5] H. Luecke, B. Schobert, H.T. Richter, J.P. Cartailler, and J.K.  Lanyi, Structure of bacteriorhodopsin at 1.55 A resolution, J. Mol. Biol. 291, 899–911 (1999).

[6] L.  Essen, R.  Siegert, W.D.  Lehmann, and D. Oesterhelt, Lipid patches in membrane protein oligomers: crystal structure of the  bacteriorhodopsin-lipid complex, Proc. Natl. Acad. Sci. U. S. A. 95, 11673–11678 (1998).

[7] U. Haupts, J. Tittor, and D. Oesterhelt, Closing in on bacteriorhodopsin: progress in understanding the  molecule, Annu. Rev. Biophys. Biomol. Struct. 28, 367–399 (1999).

[8] J.K.  Lanyi, Molecular mechanism of ion transport in bacteriorhodopsin: insights from crystallographic, spectroscopic, kinetic, and mutational studies, J. Phys. Chem. B 104, 11441–11448 (2000).

[9] J.K.  Lanyi, Bacteriorhodopsin, Annu. Rev. Physiol. 66, 665–688 (2004).

[10] J.K.  Lanyi, Proton transfers in the  bacteriorhodopsin photocycle, Biochim. Biophys. Acta 1757, 1012–1018 (2006).

[11] A.  Kawanabe, Y.  Furutani, K.H.  Jung, and H. Kandori, Photochromism of Anabaena sensory rhodopsin, J. Am. Chem. Soc. 129, 8644–8649 (2007).

[12] H. Kandori, Ion-pumping microbial rhodopsins, Front. Mol. Biosci. 2, 52 (2015).

[13] C. Wickstrand, R. Dods, A. Royant, and R. Neutze, Bacteriorhodopsin: Would the real structural intermediates please stand up?, Biochim. Biophys. Acta Gen. Subj. 1850, 536–553 (2015).

[14]L.S.  Brown, L.  Bonet, R.  Needleman, and J.K. Lanyi, Estimated acid dissociation constants of the  Schiff base, Asp-85, and Arg-82 during the bacteriorhodopsin photocycle, Biophys. J. 65, 124–130 (1993).

[15] A. Accardi and C. Miller, Secondary active transport mediated by a  prokaryotic homologue of ClC Cl channels, Nature 427, 803–807 (2004).

[16] Y.  Guo, F.E.  Beyle, B.M.  Bold, H.C.  Watanabe, A. Koslowski, W. Thiel, P. Hegemann, M. Marazzi, and M. Elstner, Active site structure and absorption spectrum of channelrhodopsin-2 wild-type and C128T mutant, Chem. Sci. 7, 3879–3891 (2016).

[17] A.Y.  Smirnov, L.G.  Mourokh, and F.  Nori, Kinetics of proton pumping in cytochrome c oxidase, J. Chem. Phys. 130, 235105 (2009).

[18] L.  Mourokh and M.  Vittadello, Physical model of proton-pumping Q-cycle in respiratory and photosynthetic electron transport chains, Chem. Phys. 530, 110638 (2020).

[19] A.  Onufriev, A.  Smondyrev, and D.  Bashford, Proton affinity changes driving unidirectional proton transport in the bacteriorhodopsin photocycle, J. Mol. Biol. 332, 1183–1193 (2003).

[20] N.  Calimet and G.  Matthias Ullmann, The  influence of a  transmembrane pH gradient on protonation probabilities of bacteriorhodopsin: The structural basis of the back-pressure effect, J. Mol. Biol. 339, 571–589 (2004).

[21] R.B. Gennis and T.G. Ebrey, Proton pump caught in the act, Science 286, 252–253 (1999).

[22] H. Luecke, B. Schobert, H.T. Richter, J.P. Cartailler, and J.K. Lanyi, Structural changes in bacteriorhodopsin during ion transport at 2 angstrom resolution, Science 286, 255–260 (1999).

[23] E.  Nango, A.  Royant, M.  Kubo, T.  Nakane, C.  Wickstrand, T.  Kimura, T.  Tanaka, K.  Tono, C. Song, R. Tanaka, et al., A three-dimensional movie of structural changes in bacteriorhodopsin, Science 354, 1552–1557 (2016).

[24] D.  Bashford and M.  Karplus, Multiple-site titration curves of proteins: an analysis of exact and approximate methods for their calculation, J. Phys. Chem. B 95, 9556–9561 (1991).

[25] J.  Iles-Smith, N.  Lambert, and A.  Nazir, Environmental dynamics, correlations, and the  emergence of noncanonical equilibrium states in open quantum systems, Phys. Rev. A 90, 032114 (2014).

[26] J.  Iles-Smith, A.G.  Dijkstra, N.  Lambert, and A. Nazir, Energy transfer in structured and unstructured environments: Master equations beyond the Born-Markov approximations, J. Chem. Phys. 144, 044110 (2016).

[27] B.P. Kietis, M. Mačernis, J. Šulskus, and L. Valkūnas, Estimation of the permanent dipole moment of bacteriorhodopsin, Lith. J. Phys. 50, 451–462 (2010).

[28] B.F.E. Curchod and T.J. Martínez, Ab initio nonadiabatic quantum molecular dynamics, Chem. Rev. 118, 3305–3336 (2018).

[29] J.K. Yu, R. Liang, F. Liu, and T.J. Martínez, First-principles characterization of the elusive I fluorescent state and the structural evolution of retinal protonated Schiff base in bacteriorhodopsin, J. Am. Chem. Soc. 141, 18193–18203 (2019).

[30] C. Lee, S. Sekharan, and B. Mertz, Theoretical insights into the mechanism of wavelength regulation in blue-absorbing proteorhodopsin, J. Phys. Chem. B 123, 10631–10641 (2019).

[31] R.  Tripathi, H.  Forbert, and D.  Marx, Settling the  long-standing debate on the  proton storage site of the  prototype light-driven proton pump bacteriorhodopsin, J. Phys. Chem. B 123, 9598– 9608 (2019).

[32] Y.  Guo, F.E.  Wolff, I.  Schapiro, M.  Elstner, and M.  Marazzi, Different hydrogen bonding environments of the retinal protonated Schiff base control the  photoisomerization in channelrhodopsin-2, Phys. Chem. Chem. Phys. 20, 27501–27509 (2018).

[33] Y.  Zhang, P.  Xie, X.  He, and K.  Han, High-efficiency microiterative optimization in QM/MM simulations of large flexible systems, J. Chem. Theory Comput. 12, 4632–4643 (2016).

[34] S. Hayashi and I. Ohmine, Proton transfer in bacteriorhodopsin: structure, excitation, IR spectra, and potential energy surface analyses by an ab initio QM/MM method, J. Phys. Chem. B 104, 10678–10691 (2000).

[35] S.  Hayashi, E.  Tajkhorshid, and K.  Schulten, Molecular dynamics simulation of bacteriorhodopsin’s photoisomerization using ab initio forces for the  excited chromophore, Biophys. J. 85, 1440–1449 (2003).

[36] C. Punwong, J. Owens, and T.J. Martínez, Direct QM/MM excited-state dynamics of retinal protonated Schiff base in isolation and methanol solution, J. Phys. Chem. B. 119, 704–714 (2015).

[37] J.-Y.  Hasegawa, K.J.  Fujimoto, and T.  Kawatsu, A configuration interaction picture for a molecular environment using localized molecular orbitals: the excited states of retinal proteins, J. Chem. Theory Comput. 8, 4452–4461 (2012).

[38] S. Zhu, M.F. Brown, and S.E. Feller, Retinal conformation governs pKa of protonated Schiff base in rhodopsin activation, J. Am. Chem. Soc. 135, 9391–9398 (2013).

[39] N.J.A.  Coughlan, K.J.  Catani, B.D.  Adamson, U. Wille, and E.J. Bieske, Photoisomerization action spectrum of retinal protonated Schiff base in the  gas phase, J. Chem. Phys. 140, 164307 (2014).

[40] M. Huix-Rotllant, M. Filatov, S. Gozem, I. Schapiro, M.  Olivucci, and N.  Ferré, Assessment of density functional theory for describing the correlation effects on the ground and excited state potential energy surfaces of a  retinal chromophore model, J. Chem. Theory Comput. 9, 3917– 3932 (2013).

[41] A.-N. Bondar, M. Elstner, S. Suhai, J.C. Smith, and S. Fischer, Mechanism of primary proton transfer in bacteriorhodopsin, Structure 12, 1281–1288 (2004).

[42] A.-N. Bondar, S. Fischer, J.C. Smith, M. Elstner, and S. Suhai, Key role of electrostatic interactions in bacteriorhodopsin proton transfer, J.  Am. Chem. Soc. 126, 14668–14677 (2004).

[43] K.  Edman, A.  Royant, G.  Larsson, F.  Jacobson, T. Taylor, D. van der Spoel, E.M. Landau, E. Pebay-Peyroula, and R.  Neutze, Deformation of helix C in the  low temperature L-intermediate of bacteriorhodopsin, J. Biol. Chem. 279, 2147–2158 (2004).

[44] G. Mathias and D. Marx, Structures and spectral signatures of protonated water networks in bacteriorhodopsin, Proc. Natl. Acad. Sci. 104, 6980–6985 (2007).

[45] S. Gozem, P.J.M. Johnson, A. Halpin, H.L. Luk, T.  Morizumi, V.I.  Prokhorenko, O.P.  Ernst, M. Olivucci, and R.J.D. Miller, Excited-state vibronic dynamics of bacteriorhodopsin from two-dimensional electronic photon echo spectroscopy and multiconfigurational quantum chemistry, J. Phys. Chem. Lett. 11, 3889–3896 (2020).

[46] D.  Friedrich, F.N.  Brünig, A.J.  Nieuwkoop, R.R.  Netz, P.  Hegemann, and H.  Oschkinat, Collective exchange processes reveal an active site proton cage in bacteriorhodopsin, Commun. Biol. 3, 4 (2020).

[47] G.  Nass Kovacs, J.-P.  Colletier, M.L.  Grünbein, Y.  Yang, T.  Stensitzki, A.  Batyuk, S.  Carbajo, R.B. Doak, D. Ehrenberg, L. Foucar, et al., Three-dimensional view of ultrafast dynamics in photoexcited bacteriorhodopsin, Nat. Commun. 10, 3177 (2019).

[48] V.R.I. Kaila, R. Send, and D. Sundholm, The effect of protein environment on photoexcitation properties of retinal, J. Phys. Chem. B. 116, 2249– 2258 (2012).

[49] M. Macernis, J. Sulskus, C.D.P. Duffy, A.V. Ruban, and L. Valkunas, Electronic spectra of structurally deformed lutein, J. Phys. Chem. A. 116, 9843– 9853 (2012).

[50] M.J. Frisch, G.W. Trucks, H.B. Schlegel, G.E. Scuseria, M.A. Robb, J.R. Cheeseman, G. Scalmani, V.  Barone, G.A.  Petersson, H.  Nakatsuji, et  al., Gaussian 16, Revision C.01, Wallingford, CT (2016).

[51] D. Oesterhelt, The structure and mechanism of the family of retinal proteins from halophilic archaea, Curr. Opin. Struct. Biol. 8, 489–500 (1998).

[52] J.K. Lanyi, Proton translocation mechanism and energetics in the light-driven pump bacteriorhodopsin, Biochim. Biophys. Acta Bioenerg. 1183, 241–261 (1993).

[53] J.K.  Lanyi, Mechanism of ion transport across membranes. Bacteriorhodopsin as a prototype for proton pumps, J. Biol. Chem. 272, 31209–31212 (1997).

[54] S.P.  Balashov, E.S.  Imasheva, R.  Govindjee, M.  Sheves, and T.G.  Ebrey, Evidence that aspartate-85 has a higher pKa in all-trans than in 13-cisbacteriorhodopsin, Biophys. J. 71, 1973–1984 (1996).

[55] R.  Rammelsberg, G.  Huhn, M.  Lübben, and K.  Gerwert, Bacteriorhodopsin’s intramolecular proton-release pathway consists of a  hydrogen-bonded network, Biochemistry 37, 5001–5009 (1998).

[56] H. Luecke, H.T. Richter, and J.K. Lanyi, Proton transfer pathways in bacteriorhodopsin at 2.3 angstrom resolution, Science 280, 1934–1937 (1998).

[57] H.-T.  Richter, L.S.  Brown, R.  Needleman, and J.K. Lanyi, A linkage of the pKa’s of asp-85 and glu-204 forms part of the reprotonation switch of bacteriorhodopsin, Biochemistry 35, 4054–4062 (1996).

[58] M. Hatanaka, H. Kandori, and A. Maeda, Localization and orientation of functional water molecules in bacteriorhodopsin as revealed by polarized Fourier transform infrared spectroscopy, Biophys. J. 73, 1001–1006 (1997).

[59] L.A. Curtiss, K. Raghavachari, P.C. Redfern, and J.A. Pople, Assessment of Gaussian-3 and density functional theories for a larger experimental test set, J. Chem. Phys. 112, 7374–7383 (2000).

[60] M. Olivucci, T. Tran, G.A. Worth, and M.A. Robb, Unlocking the double bond in protonated Schiff bases by coherent superposition of S1 and S2, J. Phys. Chem. Lett. 12, 5639–5643 (2021).

[61] T.  Tran, G.A.  Worth, and M.A.  Robb, Control of nuclear dynamics in the  benzene cation by electronic wavepacket composition, Commun. Chem. 4, 48 (2021).

[62] H. Kandori, Y. Yamazaki, J. Sasaki, R. Needleman, J.K. Lanyi, and A. Maeda, Water-mediated proton transfer in proteins: An FTIR study of bacteriorhodopsin, J. Am. Chem. Soc. 117, 2118–2119 (1995).

[63] R.  Neutze, E.  Pebay-Peyroula, K.  Edman, A. Royant, J. Navarro, and E.M. Landau, Bacteriorhodopsin: a  high-resolution structural view of vectorial proton transport, Biochim. Biophys. Acta Biomembr. 1565, 144–167 (2002).

[64] J.K. Lanyi and B. Schobert, Structural changes in the  L photointermediate of bacteriorhodopsin, J. Mol. Biol. 365, 1379–1392 (2007).

[65] E.  Nango, A.  Royant, M.  Kubo, T.  Nakane, C.  Wickstrand, T.  Kimura, T.  Tanaka, K.  Tono, C. Song, R. Tanaka, et al., A three-dimensional movie of structural changes in bacteriorhodopsin, Science 354, 1552–1557 (2016).

[66] P. Kietis, M. Vengris, and L. Valkunas, Electrical-to-mechanical coupling in purple membranes: membrane as electrostrictive medium, Biophys. J. 80, 1631–1640 (2001).


M.Mačernis a, L. Mourokh b, L.Vinciūnas a, L.Valkūnas a,c

aVilniaus universiteto Fizikos fakulteto Cheminės fizikos institutas, Vilnius, Lietuva

bNiujorko miesto universiteto Kvinso koledžo Fizikos skyrius, Niujorkas, JAV

cFizinių ir technologijos mokslų centro Molekulinių darinių fizikos skyrius, Vilnius, Lietuva


Bakteriorodopsine (BR) vykstanti protono pernaša plačiai nagrinėjama daugeliu eksperimentiniais ir teoriniais metodais. Pastarieji tyrimai paremti ypač detaliais kristalografiniais BR duomenimis. Vis dėlto, nepaisant turimų duomenų, retinalyje vykstanti protono pernaša, inicijuota elektroninio sužadinimo, nėra visiškai aiški. Naudodami kvantinės chemijos skaičiavimus, nagrinėjome retinalio molekulę, kuri yra fiksuota 13-cis,15-syn konfigūracijos pagrindinėje būsenoje ir susidaro dėl baltyminės aplinkos įtakos. Naudodami šią fiksuotą būseną, nustatėme, kad protono pernaša į asparto rūgštį Nr.  85 (Asp-85) vyksta tarpininkaujant vandens molekulei, o analogiškoje aplinkoje taip pat dalyvaujama vėl protonuojant Šifo bazę. Aptarėme BR vykstančio protonų perdavimo modeliavimo iššūkius ir metodus, analizuodami procesą, kai molekulės elektroninio sužadinimo savybes nulemia baltyminė aplinka.