INTRODUCTION

TOP ABSTRACT INTRODUCTION METHODS RESULTS AND DISCUSSION CONCLUSIONS REFERENCES

The response of living cells to a wide range of extracellular stimuli, including hormones and neurotransmitters, is mediated by a superfamily of membrane receptors known as G-protein-coupled receptors (GPCRs). The broad spectrum of signals sensed by GPCRs makes them one of the most intriguing targets for pharmacological interventions. All GPCRs exist in equilibrium between their active and inactive forms. In the absence of an agonist, the equilibrium is on the side of the quiescent form. Binding of the agonist shifts the equilibrium toward the active form, thereby triggering the activation of the G-protein (Gether and Kobilka, 1998; Lean et al., 1980), which is bound to the cytoplasmic side of the receptor.

Rhodopsin (Rh), a photoreceptor abundantly present in the outer segment of rod cells, is the only member of GPCRs for which a high-resolution crystal structure is available (Palczewski et al., 2000; Teller et al., 2001; Okada et al., 2002). It is composed of a seven-transmembrane helical protein, opsin (~40 kDa), and of retinal, a polyene chromophore bound to the protein through a protonated Schiff base (PSB) linkage to Lys-296. The counterion for the positively charged Schiff base group is a highly conserved glutamate residue, Glu-113 (Sakmar et al., 1989), which stabilizes the PSB by increasing its pK_a (Zhukovsky and Oprian, 1989; Ebrey, 2000) and is responsible for a major part of the bathochromic shift to 500 nm of the absorption maximum of retinal in Rh.

In the dark state, retinal is in its 11-cis form and acts as an antagonist, i.e., it stabilizes the inactive conformation of Rh. Photoinduced isomerization of the chromophore to the all-trans form turns it into an agonist that triggers the transition of the receptor to the active state, which, in turn, catalyzes the exchange of GDP for GTP in the -subunit of transducin, the G-protein associated with Rh.

Retinal's photoisomerization initiates a photocycle during which a thermally driven series of interconversions between Rh's different conformational states happen. During the photocycle, the absorbed photon energy, initially stored internally in the form of a twisted structure of retinal, is used to induce conformational changes of retinal, of its binding pocket, and of the whole protein (Jang et al., 2001; Borhan et al., 2000; Pan and Mathies, 2001; Lewis et al., 1995; Han et al., 1997). Low-temperature and time-resolved UV/VIS spectroscopic measurements have shown that photobleaching of Rh involves several intermediates and most likely has more than one pathway (Kliger and Lewis, 1995; Pan and Mathies, 2001; Kandori et al., 2001, and references therein).

The first intermediate of Rh's photocycle is the strongly red-shifted (570 nm) photorhodopsin (PHOTO) and decays in ~5 ps to bathorhodopsin (BATHO) with an absorption maximum (529 nm) that still shows a bathochromic shift of 29 nm relative to Rh. BATHO is in equilibrium with the blue-shifted intermediate (BSI), which decays to lumirhodopsin (LUMI) in 150 ns. LUMI undergoes a transition to metarhodopsin-I (META-I), an intermediate that is in equilibrium with the G-protein activating form, metarhodopsin-II (META-II).

X-ray crystallographic studies of Rh (Palczewski et al., 2000; Teller et al., 2001; Okada et al., 2002) have revealed many structural details of the receptor in the dark state. All transmembrane helices are bent to different extents, except for helix III, which is located in the middle of Rh's helical bundle. The kinks, which are located near the retinal binding pocket at proline and glycine sites, may serve as molecular hinges that facilitate the conformational changes associated with the chromophore isomerization.

Although there is a consensus about the 11-cis configuration of retinal in the dark state, the orientation of the -ionone ring is still under debate (Singh et al., 2001). While several NMR studies suggest a 6s-cis conformation for the chromophore (Mollevanger et al., 1987; Smith et al., 1987; Creemers et al., 2002), a recent paper provides evidence in favor of a 6s-trans form (Gröbner et al., 2000). The refined crystal structure (Teller et al., 2001) allows for both orientations, although the 6s-cis form (72° for the dihedral angle C₅==C₆C₇==C₈) has been more favorably considered. It is noteworthy that Rh can be regenerated with either 6s-cis-locked or 6s-trans-locked retinal analogs (Imamoto et al., 1996; Ito et al., 1992). It seems that the chromophore binding pocket is large and flexible enough to accommodate a variety of different chromophores, as suggested by many experiments with retinal analogs (Jang et al., 2001; Lewis et al., 2001; Fujimoto et al., 2001; Han et al., 1997).

Salt bridges play important roles in the activation mechanism of Rh. A salt bridge between the PSB group and its counterion, Glu-113, stabilizes the connection between helices III and VII. Disruption of this salt bridge results in constitutive activation of Rh in the absence of the chromophore (Robinson et al., 1992). Other functionally important salt bridges exist in the cytoplasmic region; there is a conserved tripeptide (D/E)R(Y/W) at the interface of helix III and the second cytoplasmic loop (I2) of all GPCRs. The charged residues of this motif are engaged in two salt bridges in Rh, one between Glu-134 and Arg-135 in helix III, and another one between Arg-135 and Glu-247 linking the cytoplasmic ends of helices III and VI.

Several important conformational changes during the photocycle have been identified experimentally. The equilibrium between BATHO and BSI is independent of the protein environment, but can be influenced by chemical modification of retinal (Kliger and Lewis, 1995, and references therein), indicating that the BATHO-to-BSI transition only involves conformational changes of retinal. The proposed barrier for the formation of BSI is the steric interaction between H₈ and the C₅-methyl group of retinal (Lewis et al., 1995, 2001). A cross-link experiment (Borhan et al., 2000) has also suggested that the BATHO-to-LUMI transition is accompanied by a large motion of the -ionone ring away from Trp-265, which is in close proximity to the -ionone ring in the dark and the BATHO state. The decay of BSI to LUMI, however, can be affected significantly by the amino acid composition of the binding pocket, but not by the chromophore (Randall et al., 1991; Kliger and Lewis, 1995). Therefore, it seems that this step involves the adjustment of the protein matrix to the isomerized chromophore (Kliger and Lewis, 1995). In LUMI, the salt bridge between the PSB and Glu-113 is broken and retinal is believed to be fully relaxed (Pan and Mathies, 2001). Breaking the salt bridge, in turn, decreases the pK_a of the PSB, thereby facilitating the proton transfer to Glu-113 at a later stage (Cohen et al., 1992; Kuwata et al., 2001). In LUMI, the -ionone ring can be cross-linked to Ala-169 in helix IV (Borhan et al., 2000), implying the rotation of helix IV, because in the dark state Ala-169 is located on the opposite face of this helix.

The transition to META-I involves major protein conformational changes, such as relative movement of the helices. In META-II, which is in a temperature- and pH-dependent equilibrium with META-I, the Schiff base proton is transferred to Glu-113. It has been shown that this proton transfer is coupled to a proton uptake from the cytoplasmic medium (Kuwata et al., 2001), supposedly by Glu-134, because E134N mutants can constitutively activate transducin in the absence of 11-cis retinal (Cohen et al., 1992). The protonation of Glu-134 takes place in META-II upon binding of transducin and facilitates breaking of the salt bridge between this residue and Arg-135 (Kuwata et al., 2001). Because Arg-135 is also engaged in a salt bridge (with Glu-247) between helices III and VI, the larger flexibility of this long side chain would increase the mobility of these helices.

The activation process in Rh is associated with a shift of the relative orientation of helices III and VI. This view is supported by experiments in which cross-linking of these helices (i.e., with zinc-binding motif or SS bond) prevents receptor activation (Sheikh et al., 1996; Yu et al., 1995). These conformational changes of the core domain affect the intracellular loops I2 and I3, the key sites for G-protein recognition and activation (Spengler et al., 1993; Pin and Bockaert, 1995; Wess, 1997). It has been proposed that the C-terminus of the -subunit of transducin binds in a pocket formed by these loops (Bockaert and Pin, 1999). A more detailed picture of the activation process was provided by the study of magnetic dipolar interactions of spin-labeled residues located at the intracellular ends of helices III and VI (Farrens et al., 1996); the results indicate that the distance between these residues changed upon light activation and were interpreted as a clockwise rotation of helix VI viewed from the cytoplasmic side (Farrens et al., 1996).

In situ modeling of retinal isomerization in Rh and other retinal proteins has been the subject of many theoretical studies (reviewed in Warshel and Parson, 2001). The effect of the binding pocket on the isomerization was simulated using different modeling approaches. Many early studies (Warshel, 1976; Birge and Hubbard, 1980) used position constraints on different atoms of retinal to represent the restraining effect of the binding pocket. In some studies the binding pocket was modeled using a deformable environment (Warshel and Barboy, 1982), which could react to some extent to conformational relaxation of retinal after isomerization. The electrostatic effect of the binding pocket has been taken into account by dielectric models (Warshel, 1976) or by including an explicit negative charge into the calculations (Birge and Hubbard, 1980; Tallent et al., 1992). More recent investigations of retinal in situ isomerization have used a full atomic representation of the protein. However, due to the availability of better structural models, most of these studies have investigated the isomerization in bacteriorhodopsin (Humphrey et al., 1998; Hayashi et al., 2002). There is a very recent report (Röhrig et al., 2002) on the application of the high-resolution structure of Rh (Palczewski et al., 2000) in MD simulations. In another recent study, Choi et al. (2002) used restraints derived from NMR measurements to model META-II, the active Rh intermediate that, due to time scale limitations, cannot be reached by conventional MD simulations.

In the present study we performed MD simulations of the photoisomerization process of the retinal chromophore and of the subsequent conformational relaxation process of Rh. For this purpose we constructed a molecular system where the high-resolution crystallographic model of Rh is embedded in its natural environment, i.e., in a lipid bilayer and water under constant temperature/pressure conditions, which permits a faithful description of the conformational relaxation process of the protein. In the following we describe first the actual simulations, then present the results and their analysis, and finally summarize the main conclusions.

	METHODS

TOP ABSTRACT INTRODUCTION METHODS RESULTS AND DISCUSSION CONCLUSIONS REFERENCES

In this section we introduce the system simulated, the molecular dynamics simulations carried out, and the force field parameters used.

Modeling of rhodopsin

A molecular model was constructed using the refined crystal structure of Rh (Teller et al., 2001; Palczewski et al., 2000) (PDB entry 1HZX). We used one of the molecules (chain A) of a twin structure provided in the PDB file. The carbohydrate moieties connected to Asn-2 and Asn-15, all ions, detergents, and additive molecules were removed. The five closest crystal water molecules to the protein were kept (residue numbers 2000, 2002, 2004, 2011, 2014 in PDB entry 1HXZ), one of which is near retinal (2014). The missing amino acids 236-240 and 331-333 were inserted using the visualization program VMD (Humphrey et al., 1996) and were minimized with the molecular dynamics program NAMD (Kalé et al., 1999).

As shown in Fig. 1 B, retinal in the crystal structure is strongly distorted in the Schiff base region. The suggested conformation of the chromophore is unlikely to properly represent the dark state of Rh, because the salt bridge between Glu-113 and the Schiff base is broken, and the Schiff base proton points away from the counterion. Due to the large distortion, the minimization algorithm settles in a local minimum corresponding to a 14-cis conformer. To optimize the structure to a 14-trans one, retinal was subjected to a constrained minimization during the first few steps. After minimization, retinal adopted a more planar chain, as also shown in Fig. 1 B, similar to retinal's conformation in the first published crystal structure (Palczewski et al., 2000), and the salt bridge between the Schiff base and Glu-113 formed.

View larger version (97K): [in this window] [in a new window]   FIGURE 1   (A) Rhodopsin (yellow) with retinal in orange, embedded in membrane (gray, red). The system is solvated in water (blue) and was simulated using hexagonal periodic boundary conditions (unit cell marked by black lines). (B) Retinal with protonated Schiff base and the counterion, Glu-113. The refined crystal structure (yellow) is highly distorted between the Schiff base nitrogen and C14, and the Schiff base points to the wrong direction. After minimization we obtain a more planar retinal (blue, white), which is very similar to the one found in the first crystal structure (Palczewski et al., 2000). (C) Equilibrated rhodopsin with water molecules suggested by DOWSER in red, and the crystal water molecules in green.

All titratable groups in the protein were considered to be charged, with the exception of Glu-122 and Asp-83. Comparison of the ground-state FTIR spectra of native Rh and the E122Q, and D83N mutants, suggests that Glu-122 and Asp-83 are protonated in the native protein (Fahmy et al., 1993). In our model, Glu-122 formed a hydrogen bond with His-211, while Asp-83 hydrogen-bonded to water molecules in an adjacent cavity.

Membrane embedding

Solvation

Molecular dynamics simulations

All subsequent computations were performed for an NPT ensemble (constant pressure and temperature) with full electrostatics calculation (PME) and nonorthogonal periodic boundary conditions with a hexagonal unit cell (a = b = 78 Å, c = 95 Å). The size of the unit cell ensured a minimum distance of 10 Å between adjacent proteins in the lattice. The program NAMD (Kalé et al., 1999) was used for the MD simulations. First, the system was minimized and simulated for 200 ps with fixed protein coordinates to allow lipids to adjust around the protein, and water to diffuse into crevices of the protein and hydrophilic parts of the membrane. Then the entire system was minimized and equilibrated for 1 ns. During equilibration the structure was stable and retinal assumed a fairly planar conformation, which showed a mild twist distributed over the main chain of the polyene.

The equilibrated structure was used for the simulation of the chromophore's isomerization in Rh. Retinal was isomerized around the C₁₁==C₁₂ double bond by transiently switching the dihedral potential energy function of this bond from the ground state form to the "isomerization" one (see below), which resulted in an 11-trans retinal in ~150 fs. After isomerization, the dihedral potential was switched back to its ground state form, and a 10-ns simulation was carried out. The required time for 1 ns was ~4.5 days on 128 processors of the Cray T3E at the Pittsburgh Supercomputing Center.

To examine the molecular process of the proposed rotation of helix VI, we used steered molecular dynamics (SMD) (Isralewitz et al., 1997, 2001; Izrailev et al., 1998; Lu and Schulten, 2000). The rotation of this helix was enforced by applying harmonic constraints that induced rotation of all heavy atoms of the helix around the helix axis with a constant angular velocity of 3.6°/ps. The force constant of the harmonic constraint was set to 10 kcal/mol/Ų. A customized version of NAMD was used for this purpose.

Parameters

Here we summarize the force field parameters used to describe the protonated retinal Schiff base. Force constants of 430 and 550 kcal/mol/Ų were used for CC single bond and C==C double bond stretching, respectively. For bond angles and impropers in retinal's main chain, we used force constants of 95 kcal/mol/rad² and 100 kcal/mol/rad², respectively. The equilibrium bond lengths and angles were taken from ab initio calculations (Tajkhorshid et al., 1997; Tajkhorshid and Suhai, 1999), and the parameters for nonbonded interactions from Nina et al. (1995); for dihedral angles we used the parameters published in Tajkhorshid et al. (2000). The parameters outlined have been used previously in MD simulations of bacteriorhodopsin (Tajkhorshid et al., 2000).

The dihedral potentials in the chromophore chain are cosine functions having minima at 0° and 180°. Due to the strong electronic delocalization effect in retinal, particularly in the protonated form (Tajkhorshid et al., 1999), the barrier against the rotation of individual bonds along the chain differs from the barriers expected for pure single and double bonds. Fig. 2 B shows the rotational barrier height of the dihedral angles along retinal's backbone, reflecting a polarized electronic structure in the protonated retinal Schiff base. We used partial atomic charges for retinal calculated by ab initio quantum mechanical/molecular mechanical methods (Hayashi and Ohmine, 2000; Hayashi et al., 2001). Fig. 2 A depicts the charge distribution in protonated retinal Schiff base, indicating delocalization of the positive charge along retinal's polyene chain.

View larger version (33K): [in this window] [in a new window]   FIGURE 2   (A) Charge distribution along the polyene chain in the ground and excited state. The charges of all connected hydrogens and methyl groups are projected on the atom in the main chain (C or N). Inset: chemical structure of 11-cis retinal with its conventional atom numbering. (B) Dihedral barrier against the rotation around the bonds of the retinal polyene chain. The difference for single and double bonds increases with the distance from the Schiff base. Inset: dihedral potential energy of the C11==C12 bond in the ground state and the "isomerization" state.

After photoabsorption, the chromophore isomerizes around the C₁₁==C₁₂ bond, changing its conformation from 11-cis to all-trans. In the ground state the applied potential function separates these two isomers by a barrier of 35 kcal/mol. To induce the isomerization, we transiently (for 200 fs) switched the C₁₁==C₁₂ dihedral potential to one with a single minimum at 180° (inset in Fig. 2 B). In the 11-cis conformation (0°) the new potential has a maximum that is ~42 kcal/mol higher than the corresponding minimum in the ground state potential. Thus, by switching to the "isomerization" potential we instantly add this amount of energy to the system. The potential shape is described in Hayashi et al. (2002). From the seven cosine Fourier components therein we used only the three largest ones, which are sufficient to reproduce the general features of the potential energy curve, such as the plateau region observed in ab initio molecular orbital calculations (Garavelli et al., 1997, 1998; González-Luque et al., 2000). This description maps the transition from an excited-state potential through a conical intersection back to the ground state onto a one-dimensional potential energy, with the dihedral angle as the only variable.

One-dimensional simulations of isomerizations have been done for a long time with Rh and bacteriorhodopsin (Warshel, 1976; Birge and Hubbard, 1980; Birge et al., 1987) and were often focused on understanding the processes in the chromophore. Here we mainly focus on the effects of the isomerization on the protein environment. A better description would include the change of the bond length of the isomerizing bond (Garavelli et al., 1997, 1998; González-Luque et al., 2000; Warshel et al., 1991; Warshel and Chu, 2001), but it is unlikely that our simplified description of the photoisomerization dynamics causes unfavorable discrepancies on a nanosecond time scale.

Along with the switching to the "isomerization" potential for 200 fs, we also changed the charge distribution of retinal to excited-state values (Hayashi and Ohmine, 2000; Hayashi et al., 2001) (Fig. 2 A). The altered electrostatic interaction added an additional energy difference between ground and excited state of ~15 kcal/mol. Thus the total energy increase upon excitation was 57 kcal/mol, which corresponds to the absorption maximum of Rh at 500 nm.

	RESULTS AND DISCUSSION

TOP ABSTRACT INTRODUCTION METHODS RESULTS AND DISCUSSION CONCLUSIONS REFERENCES

In our MD simulations significant changes in the structure and modes of interaction of retinal and amino acids of the binding pocket, such as the relocation of the -ionone ring, were observed during the isomerization and the following 10 ns of relaxation. The isomerization completes in ~150 fs, rendering retinal in a highly twisted structure. During the first 100 ps after isomerization, retinal's structure remains the part most affected by the isomerization event, while the other parts of the binding pocket are not influenced. Within the following 600 ps, however, several dramatic structural changes of both retinal and the surrounding binding pocket occur. Examination of the 10-ns simulation after the isomerization reveals structural changes in helices and intracellular loops, which appear to be indicative of events that take place in later steps of Rh activation.

Retinal's conformation in the dark state

In the refined crystal structure of Rh, on which our model is based, retinal is highly distorted (Teller et al., 2001). Although a major part of the unfavorable distortion is corrected during the minimization and 1 ns of equilibration (see Methods), the chromophore does not assume a fully planar conformation in the dark state of Rh and maintains a significant overall twist. A ground-state twisted chromophore has also been suggested by optical spectroscopy, NMR (Lansing et al., 2002), crystallography, and MD simulations studying the ground state of other retinal proteins, such as bacteriorhodopsin (Tajkhorshid et al., 2000; Hayashi and Ohmine, 2000) and sensory rhodopsin II (Hayashi et al., 2001), and seems to be a common feature in this family of proteins.

In the dark state, the direction of retinal's backbone twist is determined by the position of the Schiff base group and its counterion, Glu-113, and the orientation of the -ionone ring. Fig. 3 depicts the twisted conformation of retinal in the dark state (phase 1). To simplify the discussion, we will use the H₇C₇N₁₆H₁₆ dihedral angle, _chain, to represent the overall twist of the retinal backbone (from the -ionone ring to the Schiff base group). In Fig. 3 (bottom) this angle is plotted against time. During phase 1, the chain twist assumes an average value of _chain = 50°, indicating a significant overall torsion.

View larger version (49K): [in this window] [in a new window]   FIGURE 3   Top: Photoinduced conformational changes of retinal and its binding pocket. Five distinct states during the conformational changes are labeled as phases 1-5: (1) the dark state; (2) from 0 ps, at which the isomerization takes place, to 100 ps; (3) 100-300 ps; (4) 300-600 ps; (5) after 600 ps. Bottom: Twist of the polyene chain in terms of the dihedral angle H7C7N16H16, chain (red) versus time for the equilibration and the first nanosecond of relaxation; nonplanarity of ring and chain measured as the dihedral angle C5==C6C7==C8 (black).

In the ground state, the largest contribution of an individual bond to the overall torsion of retinal's backbone is from the C₁₁==C₁₂ double bond ( = 13°). Due to the steric repulsion between the C₁₃-methyl and H10 groups in the 11-cis form, the C₁₁==C₁₂ double bond deviates significantly from planarity. A pre-twisted dihedral angle, in turn, can create a strong preference for retinal to isomerize always in the same direction, although the potential energy function used for the induction of isomerization is impartial with regard to the direction of isomerization. Results of a recent theoretical study in bacteriorhodopsin clearly suggest that a pre-twisted chromophore can lead to a unidirectional isomerization of retinal after photoabsorption (Hayashi et al., 2002).

Isomerization dynamics

Fig. 4 shows the dynamics of retinal after photoisomerization at t = 0, as reflected in the time evolution of its dihedral angles and various contributions to its energy. Characterized by the change of the C₁₁==C₁₂ dihedral angle shown in Fig. 4 B, the isomerization takes place within ~150 fs, which is in good agreement with experimental observations (Schoenlein et al., 1991). A similar time scale has also been obtained in one of the early simulations of retinal isomerization reported by Warshel (1976), although the isomerization dynamics, described according to a bicycle-pedal model, is different from the present study. The final dihedral angle averages around 150°. As shown in Fig. 4 C, the total energy of retinal, defined as the sum of the conformational energies and the nonbonded interactions within retinal and between retinal and its environment, exhibits a drop of ~28 kcal/mol after the isomerization and shows an oscillation with a period of 200 fs. In Fig. 4 A, the evolution of the dihedral energies is shown for every dihedral angle in the polyene chain during the first 500 fs after switching to the "isomerization" potential. The high peak for the C₁₁==C₁₂ bond represents the maximum 42 kcal/mol of the isomerization potential.

View larger version (52K): [in this window] [in a new window]   FIGURE 4   (A) Twist along retinal backbone. The propagation of twist along the retinal backbone (atoms C6 through N16) is shown in terms of dihedral energies for the first 500 fs after switching to the "isomerization" potential. After 200 fs the potential is switched back to the ground state. (B) Dihedral angle of the C11==C12 bond versus time during the same time period. (C) Dihedral energy of retinal during the first picosecond of isomerization and relaxation.

Throughout the isomerization and subsequent relaxation, the energy is distributed over several bonds with the highest contributions from the C₉==C₁₀ and C₁₁==C₁₂ double bonds, while the single bond between them contributes very little. The reasons are steric restrictions for the retinal conformation in the binding pocket, such as the interaction of C₉-methyl with Tyr-268 and the orientation of the bonds with respect to the chain twist. The only other dihedral angle that has a noteworthy energy at the beginning of the isomerization is C₅==C₆C₇==C₈, and this is due to the nonplanarity of ring and chain. In Fig. 3 (bottom), one can see that the C₅==C₆C₇==C₈ dihedral angle decreases from 50° to 30° after isomerization, thereby accommodating some of the torsion introduced into the chain by the isomerization.

Immediately after isomerization, retinal adopts a much higher degree of torsion in its backbone (_chain = 100°) as shown in Fig. 3 (bottom). The highly twisted chromophore relaxes during the next 600 ps in several steps to a significantly less distorted structure (_chain = 150°). We note that the sign of this angle does not change after isomerization. Considering the conversion of 11-cis retinal to its all-trans form, the sign conservation of _chain indicates that the handedness of the backbone twist changes during the isomerization.

During the first 100 ps after isomerization (phase 2 in Fig. 3), the amino acids surrounding retinal remain unaltered. To still fit into the binding pocket, retinal's polyene chain kinks, i.e., retains its overall L-shaped structure. The kink is achieved by rotation around the C₉==C₁₀ and C₁₁==C₁₂ bonds, which also account for the high energy values in Fig. 4 A. A similar twisted structure was found in an earlier theoretical study investigating the isomerization of retinal in a constraining, but deformable, cavity (Warshel and Barboy, 1982). The authors proposed an intermediate with 150° and 140° for the C₉==C₁₀ and C₁₁==C₁₂ dihedrals, respectively, while we obtain values of 144° and 151°. They estimate a strain energy of ~25 kcal/mol, while we get a value of ~17 kcal/mol (sum of vdW and bonded energies). The tension related to this high conformational energy is compensated by the face-to-face van der Waals (vdW) interaction of the -ionone ring with Trp-265, as seen below. The effect of the constraints imposed by the binding pocket on the shape of the resonance Raman spectra of bathorhodopsin was analyzed by Warshel and Barboy (1982), where certain peaks in the spectra were related to a strained conformation of all-trans retinal. Phase 2 of our simulation corresponds to these observations.

In the next section we will see how the energy stored in this conformation drives further relaxation processes that involve significant changes of the binding pocket and during 600 ps leads in several steps to a much less twisted retinal (_chain = 150°).

Binding pocket relaxation

Before the isomerization, Trp-265 is positioned among the -ionone ring, the C₁₃-methyl group of retinal, and the side chain of Phe-261 (phase 1 in Fig. 3). The aromatic ring of Trp-265 is parallel to the -ionone ring. The face-to-face contact between the -ionone moiety and Trp-265 hinders the rotation of either rings during the equilibration, in accord with a recent NMR study (Creemers et al., 2002). Although there seems to be sufficient space available for fluctuation of this amino acid (Fig. 3), the strong attractive vdW interaction between the two rings reduces the chance of Trp-265 to sample this free volume. Trp-265 also provides an anchor for helix VI, which has the least interaction with other helices through hydrogen bonds, namely only between Tyr-268 and Glu-181 in extracellular loop E2 and through a weak hydrogen bond between Met-257 and the water cluster Wat1a-e mentioned in the Methods section. This relative lack of helix-helix interaction leaves helix VI free to rotate in a later relaxation step, as we will discuss further below.

Upon isomerization, however, C₁₃-methyl is relocated, significantly decreasing its contact with Trp-265, which removes one of the constraints for the position of this tryptophan. As shown in Fig. 3 (phase 3), a significant conformational change takes place in the binding pocket after 100 ps. The high tension of the kinked and twisted retinal is released by straightening the bend of the retinal chain, thereby moving the -ionone ring away from Trp-265 and closer to Met-207 (not shown in Fig. 3) in helix V, pushing it ~1 Å aside. Meanwhile, the chain twist also relaxes slightly. These changes in retinal are accompanied by a rotation of Tyr-268 to make space for the C₉-methyl group in the new conformation. This process is described in further detail below.

Trp-265 no longer impedes the rotation of the -ionone ring, and after 300 ps the ionone ring flips over, leading to a distorted 6s-trans chromophore (phase 4 in Fig. 3). The C₅==C₆C₇==C₈ dihedral angle, illustrated in Fig. 3 (bottom), shows this transition clearly as it changes from 40° to 150°. Now H₈ sits between the two C₁-methyl groups, which are both in contact with one face of Trp-265. At the same time, retinal adopts again a higher overall torsion along the main chain, _chain = 100° (Fig. 3). Finally, through a concerted rotation of ring and chain after 600 ps, retinal reaches a ring-chain coplanar, 6s-trans conformation with only a slight torsion of the chain (phase 5 in Fig. 3).

The space taken by the -ionone ring after isomerization is occupied in the dark state by the side chains of Met-207, Cys-167, and Thr-118, and a water molecule. During the relaxation of the isomerized retinal, these residues are slightly displaced to the side, to make room for the -ionone ring.

The distance between C₃ of the -ionone ring and Ala-169 in helix IV decreases ~3-4 Å during the relaxation processes. In the dark state, the -ionone ring lies further to the cytoplasmic side of helix IV and moves 3-4 Å toward the extracellular side after isomerization and relaxation. The decrease of the distance of Ala-169 and the -ionone ring was also observed in another MD study (Röhrig et al., 2002) and is in line with experimental evidence that Ala-169 is the cross-linking partner of the -ionone ring in LUMI and later intermediates, whereas Trp-265 is the preferred partner in the dark state (Borhan et al., 2000). The relocation of the -ionone ring is toward helix IV and brings it already closer to Ala-169. However, further conformational changes, e.g., rotation of helix IV around its axis, would be necessary for a successful cross-linking between the two groups, because Ala-169 is located on the opposite side of the helix. Such large conformational changes can only happen on a much longer time scale than covered in our simulations.

As mentioned above, the loss of the face-to-face interaction between the -ionone ring and Trp-265 in phase 3 is accompanied by a large rotation of Tyr-268 by 70°. Fig. 5 A clearly shows the orientational change of Tyr-268 in the binding pocket during the transition from phase 2 to phase 3. This residue is located in helix III just one pitch toward the extracellular side from Trp-265. Before the transition, both the C₉- and the C₁₃-methyl groups interact with the same face of Tyr-268, whereas after the rotation they sandwich the hydroxyl group. The rotation is induced by steric interaction with the methyl groups and by a slight motion of the polyene chain away from helix VI that accompanies the separation of the -ionone ring and Trp-265. Upon isomerization, the vdW interaction energy between Tyr-268 and retinal rises by ~4 kcal/mol, due to the relocation of the methyl groups. Tyr-268 cannot simply retreat from this steric interaction because it is hydrogen-bonded to Glu-181, which restricts the rotation of Tyr-268.

View larger version (51K): [in this window] [in a new window]   FIGURE 5   Structural changes induced by retinal isomerization. (A) The 100 ps transition in the relaxation of retinal is accompanied by a rotation of Tyr-268 to make space for the C9-methyl group in the new conformation. Left and right are the structures of the binding pocket before and after the transition, respectively. (B) About 550 ps after isomerization a change occurs in the hydrogen-bonding network of the binding pocket, including the Glu-113 counterion of the retinal Schiff base. Before the transition (left), Glu-113 in helix III is bound to Cys-187 in the -sheet segment of the extracellular loop E2, and Thr-94 in helix II hydrogen bonds to Cys-185 in the same loop. After the transition (right) Glu-113 abandons Cys-187 and bonds to the backbone and side chain of Thr-94 which, in turn, also loses its bond with the loop. The new interaction of Thr-94 interferes with the intrahelical hydrogen bonding of helix II and increases the kink. The hydrogen bonds Gly-89Thr-93 and Gly-90Thr-94 are lost and replaced by the one between Gly-90 and Thr-93.

During phase 4 the interaction energy drops again because the position of the methyl groups again resembles the one in phase 2, causing Tyr-268 to rotate back. This shows a close coupling between the orientations of Tyr-268 and the polyene chain. Finally, in phase 5, the vdW interaction energy reaches a value that is ~4 kcal/mol higher than before the isomerization, and again Tyr-268 rotates to assume the same orientation as in phase 3. At the same time we measured stabilization in the electrostatic interaction of ~5 kcal/mol because of the attraction between the negative partial charge of the hydroxyl oxygen and the delocalized positive charge of the protonated Schiff base group. We will analyze later the energetics of the conformational changes.

Another significant structural change arises around Glu-113, i.e., the counterion of the protonated Schiff base, at the end of phase 4, 556 ps after isomerization. Fig. 5 B shows the hydrogen-bond network around Glu-113. Before the structural change, helix III is connected to the -sheet on the extracellular side of the protein through a hydrogen bond between Glu-113 and Cys-187 (Fig. 5 B, left). At the same time a hydrogen bond between Thr-94 and Cys-185 connects helix II to the -sheet. In the transition these hydrogen bonds are replaced by bonds between Glu-113 and both the backbone amino group and the hydroxyl group of Thr-94 (Fig. 5 B, right). The intrahelical hydrogen bonds in helix II also undergo a reordering at the glycine kink of Gly-89/Gly-90. Before the reordering, Gly-89 and Gly-90 were hydrogen-bonded to Thr-93 and Thr-94, respectively. After the hydrogen-bond rearrangement, however, Thr-94 and Thr-93 form hydrogen bonds with Glu-113 and Gly-90, respectively, leaving the carbonyl group of Gly-89 with no hydrogen-bond partner. These rearranged intrahelical hydrogen bonds should increase the flexibility of helix II in this region and might thus provide a hinge for later conformational changes. Our findings are in agreement with magnetic dipolar interaction measurements suggesting an outward movement of the cytoplasmic part of helix II (Altenbach et al., 2001).

As previously discussed, after isomerization retinal is under mechanical stress and relaxes within 600 ps. To clarify the energetics of the process, we analyzed retinal's conformational energy and its interaction energy with the surroundings. Fig. 6 shows the changes in the various energy contributions during the 10 ns of simulation. First the dihedral energy rises because of the increased twist in the chain, but when reaching phase 5, retinal attains an almost planar conformation and the dihedral energy drops to a value 3 kcal/mol lower than for the initial configuration. Another 3 kcal/mol is contributed by the relaxation of the bond angles, while the improper energy does not change significantly. The altered interactions of retinal and its environment increase the vdW energy by ~6 kcal/mol, mostly due to loss of the face-to-face interaction between the -ionone ring and Trp-265. The electrostatic energy also rises by ~6 kcal/mol. This analysis clearly indicates that the energy initially stored in the conformational distortion of retinal is transformed into a change of the nonbonding interaction of retinal with its environment.

View larger version (30K): [in this window] [in a new window]   FIGURE 6   Bonded and nonbonded energies of retinal during the entire simulation. All energies are measured relative to initial values. Isomerization leads to a transient increase of the dihedral energy while the angle energy drops immediately after isomerization. Both energies stabilize at a value of ~3 kcal/mol. The improper energy is nearly unaffected, but van der Waals and electrostatic energies rise after isomerization and settle at values of ~+6 kcal/mol. After 10 ns of relaxation the total energy content of retinal is +7 kcal/mol.

Protein conformational change

On the nanosecond time scale some motions of the -helices and a significant movement of the intracellular loop I3 (see below) can be observed. Other major conformational changes were not clearly seen within the 10-ns simulation, although those changes are also expected to take place in the signaling process. Nevertheless, one can discern that the helices exhibit different patterns of mobility that might be related to the conformational changes on a longer time scale. The root-mean-square displacement (RMSD) of the transmembrane helices is depicted in Fig. 7. As a reference, a structure taken after 500 ps of equilibration was used. The RMSD values were evaluated for only the backbone atoms of the helices. The RMSD for each helix was computed while using the other six transmembrane helices for structural alignment. Helices VI, V, and I display the highest deviations from the reference structure. In another molecular dynamics study of Rh's photocycle, the authors report the largest RMSDs for helices IV, V, and VI (Röhrig et al., 2002).

View larger version (69K): [in this window] [in a new window]   FIGURE 7   Left: RMSD of the transmembrane helices. The reference structure for the alignment was taken after 500 ps of equilibration, using the backbone of all transmembrane helices except the one the RMSD is computed for, respectively. The largest drift can be seen for helices I, II, V, and VI, while helices IV and VII are stable and helix III shows only minor movement. Right: Superposition of the averaged structure before and 10 ns after isomerization. Coordinates were averaged over 250 ps. The color of each atom corresponds to the distance between the two positions. Blue indicates small and red large distance. The cytoplasmic region is very mobile, especially the loop connecting helices V and VI and the C-terminus. Together with loop I3, the cytoplasmic parts of helices V and VI are displaced. The extracellular half of helix I drifts ~2 Å, which is represented in an RMSD of 1.3 Å for the entire helix. The kink region displays increased mobility in helices II (behind I), VI (center, background), VI, and VII. Note that the RMSD values of the left and right figures are not directly related because of the different alignment procedures.

The right side of Fig. 7 shows the superposition of representative structures from just before and 10 ns after isomerization. The color of each segment corresponds to the extent of its displacement from the reference structure, blue being small and red large. The largest differences are in the cytoplasmic region of helix VI, most likely due to the migration of loop I3 that connects it to helix V, as described later in this section. It is hard to localize the molecular origin of the RMSD change for helix V, but the hump in the RMSD roughly synchronizes with the motion of loop I3, which suggests a connection. The extracellular end of helix I shows slight movement. Behind helix I, one can see a large movement in helix II near the Gly-89/Gly-90 hinge. Results of site-directed spin-labeling experiments were interpreted as an indication of movements of the cytoplasmic parts of helices II, VI, and VII (Altenbach et al., 2001), facilitated by the flexibility of proline and glycine hinge regions.

After ~7 ns the intracellular loop I3 starts to fold inward. This affects the structure of helices V and VI and changes the structure of the salt bridges between Glu-247 and Arg-135. Fig. 8 A compares averaged structures from before and 8 ns after the isomerization. Later in the simulation, loop I3 partially moves outward again. These movements cannot be connected to the processes in the binding pocket, but at least demonstrate the flexibility of the cytoplasmic region, which is important for G-protein binding (Arimoto et al., 2001; Bockaert and Pin, 1999; Choi et al., 2002, and references therein).

View larger version (38K): [in this window] [in a new window]   FIGURE 8   Structure of the cytoplasmic region. (A) Comparison of averaged structures of the cytoplasmic region before (yellow) and 8 ns after (blue) isomerization. After ~7.5 ns the intracellular loop I3 started to fold inward. This affected the structure of helices V and VI and altered the salt bridges between Glu-247 and Arg-135. (B) Highly conserved salt bridge Glu-134-Arg-135 and the salt bridge Arg-135-Glu-247 connecting helix III and helix VI. A typical hydrogen-bond network around Glu-134 is shown. (C) Water molecules attacking the Glu-134-Arg-135 salt bridge.

Fig. 8 B displays two important salt bridges, Glu-134-Arg-135 and Arg-135-Glu-247, in the cytoplasmic region. Glu-134 is believed to become protonated in META-II (Arnis et al., 1994; Fahmy et al., 2000), and Arg-135-Glu-247 connects Helices III and VI. According to a recent study (Choi et al., 2002), the salt bridge between Arg-135 and Glu-247 is broken in META-II. Glu-134 and Arg-135 are surrounded by hydrophobic residues toward the extracellular side, while their cytoplasmic side faces the bulk water. Water molecules can diffuse into the hydrophobic region and interact with Glu-134. Although those salt bridges were stable throughout the simulation, the one between Glu-134 and Arg-135 is challenged for a moment by several water molecules that form a complex with one of the carboxyl oxygens of Glu-134, as shown in Fig. 8 C. This happened only once during the 10-ns simulation. The weakening of the salt bridge might be related to the movements of helix VI and loop I3, since the average distance between Glu-247 and Glu-134 decreased by ~1 Å (Fig. 8 A) due to those conformational changes, leading to stronger repulsion between the two charged groups. Furthermore, the infolding loop I3 brings another charged side group, Glu-239, into the vicinity of the salt bridge that alters the electrostatic environment in this region.

Protonation of Glu-134 and the proton transfer from the Schiff base to Glu-113 are coupled (Kuwata et al., 2001). The functional importance of this salt bridge was also demonstrated by experiments in which E134N mutants could constitutively activate transducin in the absence of 11-cis retinal (Cohen et al., 1992). If the salt bridge disrupts upon proton uptake, the long side chain of Arg-135 could allow for an increased flexibility of helices III and VI.

Rotating helix VI

Several experiments have revealed that a movement of helix VI relative to the other helices plays a key role in Rh activation (Farrens et al., 1996; Altenbach et al., 2001; Sheikh et al., 1996). The movement has been mainly interpreted as a rotation of the helix. The present MD simulations suggest that the rotation is linked to the substantial shift of the -ionone ring after retinal's isomerization. As seen in Fig. 9 A, in the dark state (phase 1) the -ionone ring is blocking the sidewise movement of Trp-265, which is necessary for the rotation of helix VI, whereas the straightening of retinal in phase 5 renders the -ionone ring far enough away from Trp-265 to remove this steric obstacle.

View larger version (28K): [in this window] [in a new window]   FIGURE 9   Steered molecular dynamics rotation of helix VI. (A) Helix VI rotated with constant angular velocity around the axis shown in green. The rotation was enforced counterclockwise for the equilibrated structure just before the isomerization (retinal in gold) where the -ionone ring interacts closely with Trp-265 and obstructs its way, and with a structure 5 ns after isomerization (retinal in yellow) with the ring moved away from Trp-265. (B) Torque applied to rotate helix VI with a constant angular velocity of 3.6°/ps plotted against the angle before (solid red line) and 5 ns after (dashed blue line) retinal isomerization. Only the atoms of residues 261-268 were included in the calculation of torque. The maximum torque needed is about the same in both systems, but in the case of dark state Rh the barrier is much wider and the high torque is exerted over a longer time.

The suggested rotation of helix VI cannot be observed within the time scale of equilibrium simulations because it takes place on a much longer time scale. To test putative long-time events, SMD simulations have been widely used (Isralewitz et al., 1997, 2001; Izrailev et al., 1998; Lu and Schulten, 2000; Jensen et al., 2002). We used the SMD technique to study the rotation of helix VI. For this purpose, a torque was applied to the heavy atoms of helix VI (see the Methods section). The rotational axis was defined as the one connecting the backbone center of mass of residues 265-277 to that of residues 253-268, and is displayed as a green arrow in Fig. 9 A. By using the SMD technique, we determined at each time step of the SMD simulation the torque needed to induce a constant angular velocity (0.36°/ps) rotation of helix VI in phase 1 (just before the isomerization) and in phase 5 (after 5 ns of post-isomerization relaxation). Fig. 9 B shows the torque profile along the rotational angle computed with forces acting only on residues 261-268. In both phase 1 and phase 5, the torques rise rapidly and reach their maxima within 20°. These increases of the torques originate from interactions of helix VI with the surroundings. The torque in phase 5 drops quickly after 20°, whereas that in phase 1 remains at the maximum value until 30°. Because the rotation is hindered by a steric interaction of Trp-265 in the 20-30° region in phase 1, a larger torque is needed to overcome a potential barrier due to the steric interaction. However, the -ionone ring positions away from Trp-265 in phase 5 (Fig. 9 A), enabling Trp-265 to pass through this region with lower torque applied.

	CONCLUSIONS

TOP ABSTRACT INTRODUCTION METHODS RESULTS AND DISCUSSION CONCLUSIONS REFERENCES

In the present study we performed a 10-ns MD simulation to examine the conformational relaxation process of Rh after the cis-trans isomerization of the retinal chromophore. Even though it is difficult to assign intermediate states observed during the relaxation in our simulation to the spectroscopically determined intermediates due to lack of statistics, our simulations are in accordance with several important experimental observations.

The twisted retinal structure in the dark state is preserved throughout equilibration. Isomerization first leads to a strongly twisted and bent retinal maintained by the close contact between the -ionone ring and Trp-265. The loss of this interaction after 100 ps allows retinal to relax to a much less distorted conformation in several steps within the next 500 ps. These steps involve the relocation of the -ionone ring away from Trp-265 and its rotation from a strained 6s-cis geometry to a coplanar 6s-trans conformation. The relocation of the -ionone ring decreases its distance to Ala-169 by 3-4 Å. This observation is in line with the reported cross-link between this residue and the -ionone ring in LUMI (Borhan et al., 2000), although further conformational changes, e.g., the rotation of helix IV in later steps, is necessary for the required alignment of these groups. In this regard the rotation of helix IV can be studied by SMD, in a similar manner used for the rotation of helix VI in this paper.

The results of Raman microchip spectra suggest that the BATHO-to-LUMI transition involves a relaxation of the distorted chromophore to an unstrained conformation (Pan and Mathies, 2001). The formation of LUMI, which takes place in 150 ns, is believed to drive the Schiff base group out of the Glu-113 region (Ganter et al., 1988; Pan and Mathies, 2001). Within a 10-ns relaxation in the present simulation this could not be observed, as the distance of the Schiff base and counterion is constant. However, energy analysis showed that there is still a significant amount of energy stored in the system after 10 ns, which could facilitate the later separation. The present simulation revealed significant conformational changes in the binding pocket and cytoplasmic region. The displacement of retinal's C₉-methyl upon the isomerization is the driving force for a conformational change of Tyr-268 in the binding pocket. In the Schiff base region we observed a reordering of the hydrogen-bond network including helices II and III. The cytoplasmic loop I3 displays large mobility, which affects an interhelical salt bridge close to the G-protein binding site.

Gross conformational changes of Rh, including overall movements of the helices required for the G-protein binding and activation, are not expected within 10 ns. Nevertheless, the observed conformational changes in the binding pocket and the cytoplasmic region alter the interaction of the helices and, therefore, signify the advent of the gross conformational changes. In fact, RMSD analysis shows different patterns of mobility for the helices, especially high mobility of helix VI. In the dark state Trp-265 provides an anchor for helix VI through the interaction with the -ionone ring and the C₁₃-methyl group. During the relaxation process this anchor becomes loose and the mobility of Trp-265 increases significantly. With the -ionone part of retinal moved out of the way like an opened turnstile, the rotation of helix VI (Farrens et al., 1996) should be alleviated.

Resolving the physical mechanism of the light activation in Rh is obviously one of the most appealing challenges in quantum biology today, now that a structure of this protein is finally available. The activation involves processes stretching over many time scales, starting with a femtosecond excitation process that triggers a 100-fs isomerization, which links into a cascading of retinal's internal energy into energy of the binding site. The energy is then utilized to prepare the protein for large-scale conformational transitions needed for Rh to interact with transducin, which eventually induces a neural signal. Our study shows that on the nanosecond time scale of molecular dynamics simulations one can catch a glimpse of the processes that activate Rh. The next round of investigations should address on the one side in how far the simple model used to describe the photoisomerization in our description is appropriate. For this purpose the dynamics of in situ retinal on its excited-state potential energy surface and the crossing to the ground state surface should be described in complete detail. On the other side, simulations should be extended to the 100 ns and microsecond range to identify the transitions that lead to interaction of Rh with transducin. The structure of the latter is known and, hence, docking between Rh and transducin can be simulated. A recent study (Choi et al., 2002) combining modeling techniques and experimental information has proposed a model for the activated form of Rh, META-II. Further MD simulations in the microsecond range can verify the model. Such simulations can also characterize the activation pathway in detail. If the time scale of microseconds cannot yet be reached, future studies might first use steered molecular dynamics to test whether suggested conformational changes in Rh are feasible. This has already been attempted in this paper for the rotation of helix VI and should be repeated using weaker external torques, avoiding strong friction. Such calculations can determine the steps involved in the activation of Rh and verify the structure of the active form, META-II.