Next Article in Journal
Analysis and Test of Internal Blowing Anti-Tangle Bag-Breaking Device for Domestic Waste
Previous Article in Journal
Optimization of the Automated Production Process Using Software Simulation Tools
 
 
Font Type:
Arial Georgia Verdana
Font Size:
Aa Aa Aa
Line Spacing:
Column Width:
Background:
Article

The Mechanism of Channel Opening of Anion Channelrhodopsin GtACR1: A Molecular Dynamics Simulation

1
Chongqing Key Laboratory of Big Data for Bio Intelligence, Chongqing University of Posts and Telecommunications, Chongqing 400065, China
2
Department of Chemistry and Physical Sciences, Nicholls State University, P.O. Box 2022, Thibodaux, LA 70310, USA
*
Author to whom correspondence should be addressed.
Processes 2023, 11(2), 510; https://doi.org/10.3390/pr11020510
Submission received: 22 December 2022 / Revised: 2 February 2023 / Accepted: 4 February 2023 / Published: 8 February 2023

Abstract

:
Guillardia theta anion channelrhodopsin 1 (GtACR1) is a widely used inhibitor of optogenetics with unique conductance mechanisms and photochemistry. However, the molecular mechanism of light-gated anion conduction is poorly understood without a crystal structure for the intermediate state. In this study, we built the dark-state model based on the crystal structure of retinal and isomerized the model by twisting the C12-C13=C14-C15 dihedral step by step using molecular dynamics simulation. The conformational changes revealed the all-trans to 13-cis photoisomerization of the retinal chromophore cannot open the channel. There is no water influx, and a pre-opened K-like intermediate after photoisomerization of retinal is formed. During the opening of the ion channel, proton transfer occurs between E68 and D234. Steered molecular dynamics (SMD) and umbrella sampling indicated that the E68 and D234 were the key residues for chloride-ion conducting. We propose a revised channel opening pathway model of GtACR1 after analyzing (de)protonation of E68 and D234. Reprotonation of D234 will result in two different early L intermediates, named L1-like and L 1 ‘-like, which correspond to the L 1 and L 1 ‘ intermediates reported in a recent study. Simulation results showed that L 1 -like may convert by parallel paths into L 1 ‘-like and L 2 -like states. This model provides conformational details for the intermediate as well.

1. Introduction

Controlling and adjusting neurophysiological activities in the cell by using light-sensitive proteins are among the significant innovations in neuroscience in recent years [1,2]. Even though channelrhodopsins (ChRs) [3,4,5] have modified neuroscience research by depolarizing membranes to direct light to activate neuro cells [6], investigations on the photoinhibition of neuronal action potentials have been limited due to the lack of effective membrane hyperpolarization tools. Govorunova et al. [7] discovered the first natural anion channel rhodopsin (ACR) from the cryptophyte alga Guillardia theta (GtACR1 and GtACR2), which is the most used optogenetic tool for neuron suppression. Since the discovery of ACRs, the gating mechanisms have been studied mainly by using absorption spectroscopy and photocurrent analysis of wildtype (WT) ACRs and their mutants [8,9,10,11,12,13,14,15,16,17,18]. Kim et al. [19] and Li et al. [20] resolved the tetramer and the dimer crystal structure of GtACR1, respectively. The high-resolution crystal structure provides atomic insights into further research on the gating mechanisms. However, accurate information on ion channel structure of ACRs is still unavailable.
The GtACR1 monomer is composed of seven transmembrane helices (TM1-7) [20]. The anion channel is surrounded by low-polarity or nonpolar aliphatic amino acids located at TM1-3 and 7 (Figure 1). The photoisomerization of retinal from all-trans to 13-cis involves a series of spectral intermediates, labeled K, L, M, and N/O. As the entrance of the ion channel, extracellular vestibule 1 (EV1) links to extracellular constriction site 1 (ECS1, including K33, I27, Y81, A75), central constriction site 1 (CCS1, including E68, N239, Q46, S43), and intracellular constriction site 3 (ICS3, including L108, L245, A61, P58, E60) successively, and a full ion tunnel form. CCS1 is a major obstacle to ion conduction when the ECS and ICS are opened [19]. Li et al. [20] considered that the ion channel of GtACR1 is more likely to pass through CCS2, which is composed of L64, M105, and T101, as seen in Figure 1. The full photocycle of GtACR1 shown in Figure 1 for the all-trans to 13-cis photoisomerization of the retinal chromophore is like that for most microbial rhodopsins [13]. The protonated retinal Schiff base (RSBH + ) is covalently linked to 13-cis retinal in the K state with strained configuration, which will relax to the L states. RSBH + is deprotonated in the M state and reprotonated to form the N state, in which the retinal is still 13-cis. The retinal is twisted into the all-trans configuration to produce the O state, which decays to the primary state. It is worth noting that the ion-conducting L state precedes the deprotonation of the RSBH + . It was proposed that both E68 and D234 were protonated in GtACR1 [13,14,19]. However, according to recent study of Dreier et al. [21] and Tsujimura et al. [22,23], D234 is more likely deprotonated in the dark state. Interestingly, a small amount of M intermediate was also formed in E68Q mutants, suggesting that the presence of E68 is not a strict prerequisite for Schiff base deprotonation, and it possibly exists as an alternative proton acceptor in GtACR1 [13]. In previous research [13], the open state of GtAcR1 was assigned to the L intermediate, and it was considered that the isomerization of retinal was sufficient to allow chloride ions to pass through the ion channel. However, the L intermediate was divided into an early non-conducting state L 1 and a late conducting state L 2 in a recent time-resolved Fourier transform infrared (FTIR) spectroscopy investigation [21]. This literature suggested that E68 was protonated in the dark state and deprotonated in the early photocycle and then reprotonated in the L 1 →L 2 transition. An additional state (L 1 ‘) was proposed during the L 1 →L 2 transition, but it is not clear whether L 1 and L 1 ‘ are one mixed state or two consecutive states. Understanding the gating mechanism is of great significance for optogenetics and requires knowledge of the degrees of freedom of each atom. In this work, we established the atomic structure model of GtACR1 based on the crystal structure and carried out molecular dynamics (MD) simulation. The local conformations were studied by analyzing the changes of hydrogen bond network and water density distribution in the ion channel.

2. Materials and Methods

2.1. System Modeling

The crystal structure of the GtACR1 dimer (PDB ID: 6EDQ) [20] was used as a template to model the monomer structure of the dark state (closed) of GtACR1 (GtAcR1-trans) by deleting half of the dimer. Based on a recent article [21,22,23], the retinal Schiff base and the key residues E68 were considered to be protonated, and D234 was considered deprotonated. The pKa calculations of three independent trajectories show that the protonation state is reasonable (as seen in Table S1). Due to no available crystal structure, the model for GtACR1 after photoisomerization (GtACR1-cis) was computationally generated as follows: the C12-C13=C14-C15 torsion in retinal was rotated in 20° steps and kept frozen, as equilibration is simulated after each step, until the dihedral angle reaches approximately 0° [24]. CHARMM-GUI [25] was utilized to generate a hydrated 1-palmitoyl-2-oleoyl-sn-glycero-3-phosphocholine (POPC) lipid membrane, in which the trans and cis GtACR1 monomer models were embedded. Sodium and chloride ions were added to neutralize each system at a concentration of 0.15 mol·L 1 . The complete simulation system consisted of the GtACR1 monomer, 115 POPC molecules, and 12,696 water molecules in a periodic box of size 70.3 × 70.3 × 126.2 Å 3 with 57,283 atoms in total. The CHARMM36m [26,27], a force field was used for proteins, lipids, and ions, and the water molecules were described by the TIP3P [28] model for water. Four protonation states of E68 and D234 in the early photocycle were investigated: E68 deprotonated, D234 deprotonated, and both residues deprotonated or protonated. The states are named E68dp-D234p, E68p-D234dp, E68dp-D234dp, and E68p-D234p, respectively. According to this naming rule, GtACR1-cis can also be called E68p-D234dp. For all six abovementioned systems, 300 ns MD simulations were run to obtain thermodynamically stable conformations. We performed at least three simulations for each system with different initial velocities and obtained similar results from the trajectories.

2.2. Molecular Dynamics Simulations

All simulations were carried out using the NAMD 2.12 package [29] using NPT in addition to periodic boundary conditions. Long-range electrostatic interactions were treated with the Particle-Mesh-Ewald [30] approach with a grid spacing of less than 1 Å. The cut-off radii of long-range electrostatic and van der Waals interactions were set to be 12 Å. Covalent bonds with hydrogen atoms were constrained using the SHAKE algorithm [31] and using an integration time step of 2 fs. Initially, the steepest descent method is used to simulate 10,000 steps to minimize the energy of the system. The Nosé–Hoover chain of thermostats was used to keep the system at 310 K. The Langevin piston algorithm was used in the heating process to keep a constant pressure (1 bar) with the damping coefficient of 1 ps 1 . Then, equilibrations were performed for 300 ns under NPT conditions with periodic boundary conditions. The system stability of molecular dynamics simulation was determined by measuring the root mean square deviation (RMSD) of all non-hydrogen atoms, using the initial structure of GtACR1 as a reference (Figure S1). Protonation states of titratable residues were determined using the PDB2PQR [32].

2.3. Steered Molecular Dynamics Simulation

The transport process is typically characterized by a substrate pathway along with associated conformational changes in the protein. However, the relevant time scales for many transport processes are beyond that afforded by MD, which is typically limited to a few microseconds currently. SMD simulations were used to gain a qualitative picture of the permeation pathway for chloride in the GtACR1 channel. We assumed constant velocity to explore possible permeation pathways through GtACR1. SMD simulations [33,34,35,36] were initialized from the representative configuration extracted from cluster analysis of the equilibration. Pulling directions were utilized to explore the pathway of chloride conductance through the putative pore of GtACR1. Retinal was held in the 13-cis conformation during the simulations. The chloride ion was placed at the extracellular entrance. The chloride ion was tagged for pulling with a spring constant of 4 kcal/mol·Å 2 and pulling velocities of 0.1 Å/ns. A Hamiltonian limit of 5 kcal/mol·Å 2 was added to the amino acid residues surrounding the seven helices of the protein to avoid any distortions of the protein as a consequence of pulling. We examined the simulation trajectory in VMD [37] and correlated the force peaks in the plot with specific events.

2.4. Umbrella Sampling

As the forces of SMD imposed are usually orders of magnitude greater than the forces experienced in a living system, interpretation of the results is often limited to a qualitative description of the possible behavior of the substrate and protein. For a quantitative picture of transport, one must turn to more advanced methods. The substrate pathway can also be characterized by the substrate’s free energy along it, known as the potential of mean force (PMF), which dictates the speed of transport and the selectivity. The reaction coordinate used for the umbrella sampling [38,39] simulations was the same reaction coordinate used for the SMD simulations, with the z distance along the channel axis. The umbrella sampling simulation could find the equilibrium state in a set of umbrella sampling windows and obtain a relatively accurate PMF curve by reweighting constant biasing potential well along a reaction coordinate. The whole span of the reaction coordinate was subdivided into equally spaced windows in the z direction; each window was simulated independently with a width of 1 Å. For 1 Å windows, a force constant of 2.5 kcal/mol·Å 2 will give a good overlap between neighboring windows. The width of the Translocation colvar was 0.1 Å. For a simple US simulation, a colvar file was written for each window, with the center parameter set to the harmonic center for that window; then, each window was run independently. Finally, we used the weighted histogram analysis method (WHAM) [40] to construct the PMF from the sorted trajectory files with 6 ns sampling simulation in each window.

3. Results

3.1. Conformational Differences between GtACR1-trans and Its 13-cis Isomer

The local conformations near CCS1 of GtACR1-trans and GtACR1-cis obtained by MD simulations are shown in Figure 2A and Figure 2B, respectively. In the GtACR1-trans conformation, D234 and E68 are the only carboxylic residues in the vicinity of the nitrogen atom of RSBH + , which is consistent with crystal structure of GtACR119. In the GtACR1-cis conformation, the downward flip of RSBH + caused by photoisomerization of retinal disrupts the hydrogen bond networks between RSBH + and the side chains of D234.
In the GtACR1-trans conformation, E68 forms a hydrogen bond network with N239 and S43 and forms an ENS triad. This simulation conformation is consistent with crystal structure [19,20]. Similar homologous triads referred to as the central gate were also found in ChR2 (E90, N258, S63) [41] and C1C2 (E129, N297, S102) [42]. Unlike cation channelrhodopsins (CCRs), the ion channel in GtACR1 is not blocked by the ENS triad but passes through the seam along the CCS1 (Figure 2A,B). Obviously, GtACR1 has a different gating mechanism from CCRs. Distances between several important residues of CCS1 of GtACR1-trans and GtACR1-cis are shown in Figure 2C. The hydrogen bond network near the central gate is almost unchanged. Simulation results show that the D234-Y207 and D234-Y72 distances remained between 1~2 Å throughout the simulation run, which is indicative of hydrogen bonding and consistent with the crystal structure of GtACR1 [19,20]. After isomerization of retinal, the nitrogen atom of RSBH + turns toward the intracellular side and breaks the D234-RSBH + hydrogen bond, as manifested in the increasing distance between RSBH + and D234. In addition, all the distances between residues pairs mentioned in Figure 2C (except D234-RSBH + ) were unchanged, which suggests that the hydrogen bond network of CCS1 stays the same even if the retinal is isomerized. This indicates that the isomerization of retinal exclusively is not enough to allow chloride ions to pass through the ion channel. In GtACR1-trans, water molecules accumulate near the K238 and N239, and there are few water molecules distributed between E68 and D234 (Figure 3A). From the crystal structure of GtACR1, Kim et al. [19] proposed two possible explanations for stabilizing the positive charge on the RSBH + in the case where both E68 and D234 are protonated. One explanation is that polarized water could bind to the Schiff base, and another is that the partial negative charge of the carbonyl in the adjacent D234 may slightly stabilize the positive proton. However, recent research [21,22,23] suggested D234 could be a counterion for RSBH + while deprotonated in dark state. Our simulation results show that the D234 directly forms a hydrogen bond interaction with RSBH + , which is consistent with recent research. After photoisomerization, changes in hydrogen bond network allow water molecules to mainly accumulate around the Schiff base and enter the interspace between E68 and D234 in GtACR1-cis (Figure 3B). This suggests that the photoisomerization of retinal favors the influx of more water in the ion channel and forms a “pre-opened” state.

3.2. The Formation of Pre-Opened State (K-like) for GtACR1-cis

Water molecules play an important role in the opening of ion channels in ChR2 [43,44]. It has been reported that during the early stage of formation of ion channel, external water molecules flow into the channel to form a “pre-opened” state, which then continues to expand to create ion conduction [42]. The grid water density for GtACR1-trans and GtACR1-cis was calculated to determine the degree of ion channel opening and is displayed in Figure 3, and the local conformations were also inserted. The water distribution around retinal bound pocket in GtACR1-trans was found to be discontinuous, which is consistent with the water density in its crystal structure [19,20], suggesting that the ion channel in natural GtACR1 is closed. However, a few water molecules enter the CCS1 in the GtACR1-cis isomer. Although the limited number of water molecules is not enough to be considered as having a complete continuous distribution in the channel, it is certain that the distribution of water molecules in the GtACR1-cis isomer is significantly different from that of the wild type. The water molecules are enriched around retinal and D234, which increases theh hydration of retinal and D234. The GtACR1-cis is the first intermediate in the photocycle following the closed state. Although there is a lack of direct spectral evidence, it can be considered as the K intermediate, or at least, a K-like pre-opened state, which allows entry of some water molecules.

3.3. Binding Sites of Chloride Ions in E68p-D234p (trans) and E68p-D234p (cis)

Since both E68 and D234 are considered as counterions, we set both E68 and D234 to be protonated; the retinal were all-trans and 13-cis, respectively. To verify the effects of E68 and D234 on the electrostatic interactions at the binding sites and to further understand the process of chloride-ion permeation through the channel, SMD simulations and umbrella sampling were utilized to investigate the transmembrane transfer of the E68p-D234p(trans) and E68p-D234p(cis). Much important information, including the interactions between ion sand residues around the channel, and the suitable ion permeation path, were obtained by recording the variations in the drawing force required versus reaction coordinates during this process, as can be seen in Figure 4. Chloride ions were moved from the extracellular N-terminal to the intracellular C-terminal (coordinate of z from 16 to −14 Å) with constant velocity of 0.1 Å/ns. In the case of GtACR1-trans, the SMD simulation was performed to verify the key bonding site, although the channel was closed. Figure 4A shows a global minimum of traction force with value of −565.24 pN located at the z coordinate equal to 2.26 Å (labeled by letter a). This coordinate corresponds to the positions of D234, E68, and N239 of CCS1, which constitute the channel bottleneck in GtACR1-trans (see Figure 2A). This suggests that the chloride ion is trapped by the hydrogen bonding network among the three neutral residues of CCS1 rather than electrostatic attraction in GtACR1-trans. Subsequently, the chloride ion is pulled away from the bottleneck, which requires a huge force to overcome the constraint of the hydrogen bond network. The maximum force with a value of 399.18 pN at z = −2.99 Å (labeled b in Figure 4A) indicates that the chloride ion has to overcome a resistance of 964.42 pN through a channel bottleneck of length 5.25 Å. In contrast, in GtACR1-cis, there are two local minima and two maxima in the traction force of chloride ions during the stretching process. The two local minima occur at z = 7.94 Å with the traction force of −536.03 pN and z = −0.75 Å with a force of −452.56 pN, (labeled c and e, respectively). The former coordinate corresponds to the positions of D234 and T71 in the bottleneck of ECS1, and the later corresponds to the positions of in the bottleneck of CCS1. These minima in the traction force suggest potential binding sites for chloride ions. Similarly, there are the two local maxima (labeled d and f, respectively, in Figure 4B), which are associated with the minima. Thus, in the K-like intermediate, there are two bottlenecks with resistances of 832.49 and 744.82 pN, respectively; and their lengths are 3.78 Å and 3.24 Å, respectively. However, both the resistance and the lengths of the two bottlenecks of E68p-D234p(cis) are smaller than those of E68p-D234p(trans), which is consistent with the notion of electrostatic interactions at the binding sites of E68 and D234, and they can influence the transmembrane transport of Cl in GtACR1 obviously. To further investigate the energy changes during the chloride ion transfer along the channel, the potential of mean force (PMF) curve of Cl was calculated based on the SMD trajectory, and the Jarzynski equality was used to obtain equilibrium free energy. The PMF corresponds to the energy barrier that must be overcome to allow ion permeation. As shown in Figure 5, there are two main energy potential wells in the penetration process of Cl via umbrella sampling method. The PMF curve shows a 16.78 kJ/mol energy barrier at 6.14 Å and another at −4.5 Å (8.75kJ/mol). The coordinates found by the SMD simulation of the GtACR1-cis are in good agreement with the position of the protein binding site from the force change diagram of ions. When the Cl ion moves through the ion channel, the stabilizing force is developed at the two binding sites, including D234 and E68. Therefore, D234 and E68 are likely to be potential binding sites.

4. Discussion

Previous studies [11,13,16] suggest that the Schiff base of GtACR1 is protonated in the K and L intermediate states and is deprotonated to form the M intermediate during the photocycle. The only carboxylic residues near the Schiff base, E68 and D234, play important roles in channel gating and Schiff base’s proton transfer. In an electrophysiological experiment [13], no fast positive current was observed in the E68Q-D234N double mutant (with both E68 and D234 kept neutral). This and the UV–visible spectrum suggest that the RSBH + cannot deprotonate to form the M state. However, for the D234N single mutant, a fast positive current, with an amplitude that increases at higher pH, was observed. This indicates that the deprotonated E68 can accept a proton from RSBH + and that the channel current generated corresponds to formation of the conductive L state when D234 is neutral. As mentioned above, the M intermediate is also generated in E68Q mutants but with a low yield, which suggests neutral E68 could also be present in the photocycle. Our simulations suggest that exclusive retinal isomerization is not enough to open the ion channel. Therefore, proton transfers refer to D234 and E68 may be the key steps to the subsequent channel expansion. An important question arises: how does the deprotonation of E68 and/or protonation of D234 lead to channel opening?
To investigate the gating mechanism of E68 and D234, we performed 300 ns MD simulations to investigate the effects of the protonation states of the two residues. The water density distribution and numbers of hydrogen bonds and water molecules in channel for different protonation states of E68 and D234 were calculated, and the results are summarized in Figure 6. When the protonation states of E68 and D234 changed, the water density distributions were significantly different. However, the numbers of hydrogen bonds and water molecules in the channel changed significantly only when both residues were deprotonated and did not change much when one residue was deprotonated or when both residues remained neutral. In the K-like state (GtACR1-cis), the Schiff base is twisted to the other side away from D234 (seen in Figure 2B), disrupting the previous electrostatic equilibrium between RSBH + and the negative D234 of the wild type. The exposed negative D234 may reprotonate from solvent or protonated E68. Note that we only assume that D234 gets a proton from the solvent and cannot predict which residue is the proton donor. The negative D234 may get a proton from the protonated E68, forming an E68dp-D234p pattern (the water density distribution and local conformation are displayed in Figure 6B and Figure 7B, respectively). The E68, D234, and N239 form a stable hydrogen bond network after E68 transfers protons to D234. The hydrogen bond between E68 and D234 implies that there is an equilibrium between E68p-D234dp and E68dp-D234p states. Alternatively, when D234 gets a proton from the solvent and E68 stays protonated (descripted by E68p-D234p), the simulations showed that the water density distribution near E68 and D234 is denser than that in the GtACR1-cis state (compared Figure 3B and Figure 6A). However, the numbers of hydrogen bonds and water molecules of the double protonated state are not very different from those the GtACR1-cis state, which means that the water in the channel has better fluidity at this time. The reason may be that the hydrogen bond between E68 and N239 existing in the K-like state disappears in the E68p-D234p state, making water molecules move more smoothly. Differently from the E68p-D234dp state, there is not the stable hydrogen bond network among E68, D234, and N239 in the E68p-D234p state (Figure 7A), which suggests the constriction site composed of a hydrogen bond network in CCS1 disappeared. Interestingly, both E68 and D234 are protonated in the E68p-D234p state, suggesting that there may be other residues as counterions that accept the Schiff-base protons and evolve into M state subsequently. Figure 6A,B displays that that water distribution in both E68p-D234dp and E68p-D234p states is discontinuous. It indicates that both the E68p-D234dp and E68p-D234p states following the K-like state are not conductive states, although they may also be considered as certain forms of the L intermediate. In this case, protonated E68 is not favorable for the formation of a conductive state.
In contrast, deprotonation of E68 from the K-like or E68p-D234dp state is the alternative (descripted by E68dp-D234dp). As seen in Figure 6E, more fluid water in CCS1 indicates that the channel will continue to expand, and the E68dp-D234dp state is more conducive than the E68p-D234dp and E68p-D234p states for further opening of the channel. This inference is consistent with Sineshchekov’s description that the deprotonation of E68 facilitates the formation of the M state [13]. The decrease in the number of hydrogen bonds in the channel of E68dp-D234dp state (Figure 6D) also confirms that the channel is more permeable. In the water distribution for E68dp-D234dp, shown in Figure 6C, it seems that the water almost fills the whole tunnel. It is clear that deprotonation of both D234 and E68 facilitates ion conduction, and the E68dp-D234dp intermediate has the properties of a channel opening state.
In a recent article [21], Dreier et al. divided the L intermediate originally thought to be a conductive state into an early non-conducting L 1 /L 1 ‘ state and a late conducting L 2 state; whether L 1 and L 1 ‘ are a mixed/inhomogeneous state or two consecutive states is unclear. They concluded that E68 is protonated in the ground state and deprotonated in L 1 and reprotonated in the L 1 →L 2 transition. Steered molecular simulation and umbrella sampling indicated that the E68 and D234 were the key residues for chloride-ion conducting. Based on our simulation results, we propose a revised channel opening model, as seen in Figure 8. E68 is protonated while D234 is deprotonated in the dark state of GtACR1-trans, and remains ub this state after the fast photoisomerization of retinal to form the GtACR1-cis, the pre-opened K-like state. There are two possibilities for reprotonation of deprotonated D234: one is to obtain a proton from the solvent, and the other is to obtain a proton from E68, forming E68p-D234p and E68dp-D234p, respectively. The E68dp-D234p may be an L 1 -like intermediate, in which case, D234 transfers a proton to water, leading to an L 2 -like state. The E68p-D234p state, which may be assigned to the L 1 ‘-like category, was obtained from D234 getting proton from solvent in the K-like state, or from E68→D234 proton transfer in the L 1 -like state. The E68dp-D234dp state, which was a possible conductive functional intermediate in our simulations according to water distribution, is formed by D234 losing a proton in the L 1 -like state. The concomitant E68p-D234p state in the transformation of L 1 -like to L 2 -like, which reasonably explains the reprotonation of E68 in the L 1 →L 2 transition [21].

5. Conclusions

The channel opening process of the anion-conducting GtACR1 was investigated by molecular dynamics simulations. The results show that the isomerization of retinal led to no water influx and the formation of a pre-opened K-like intermediate. E68 and D234 each transfer a proton in the process of channel opening. SMD simulations and umbrella sampling analysis for trans and cis conformations of GtACR1 suggested that E68 and D234 are key residues of channel opening. Based on the results, we propose a revised photocycle model for GtACR1. E68 is likely to lose a proton, and D234 may obtain a proton from the solvent due to hydration, and early non-conducting L states form in the subsequent evolution of the photocycle, which are assigned to L 1 -like and L 1 ‘-like states. The L 1 -like state with deprotonated E68 and protonated D234 (E68dp-D234p) may convert in parallel paths into L 1 ‘-like and L 2 -like (E68dp-D234dp), respectively. Therefore, the conducting state in photocycle is a mixture of conducting L 2 -like and nonconducting L 1 ‘-like states. Although direct experimental evidence for this modified photocycle model is currently lacking, we hope that other experiments will follow to confirm this result.

Supplementary Materials

The following supporting information can be downloaded at: https://www.mdpi.com/article/10.3390/pr11020510/s1. Figure S1: The time evolution of root mean square deviation (RMSD) of non-hydrogen atomic positions of (A) GtACR1-trans, (B) GtACR1-cis, (C) E68dp-D234dp, (D) E68p-D234dp, and (E) E68dp-D234dp. Table S1: The pKa values of molecular dynamics clustering conformation of each representative conformation. Table S2 Atom type definition, charges and connectivity for retinal. Table S3 Related bond parameters used in the force field for describing the C12-C13=C14-C15 groups in the protein. Units are presented in the brackets. Table S4 Related angle parameters used in the force field for describing the C12-C13=C14-C15 groups in the protein. Units are presented in the brackets. Table S5 Related dihedral parameters used in the force field for describing the C12-C13=C14-C15 groups in the protein. Units are presented in the brackets. n indicates multiplicity.

Author Contributions

Conceptualization, S.Y.; methodology, S.Y. and Y.D.; software, C.L.; validation, C.L. and Q.X.; formal analysis, C.L.; investigation, C.L., C.Q. and M.J.; data curation, C.L. and Q.X.; writing—original draft preparation, C.L. and S.Y.; writing—review and editing, G.V.L.; visualization, C.L. and Q.X.; supervision, S.Y. All authors have read and agreed to the published version of the manuscript.

Funding

This research received no external funding.

Data Availability Statement

Not applicable.

Acknowledgments

We thank Shuangyan Zhou for her help in molecular dynamics simulation.

Conflicts of Interest

The authors declare no conflict of interest.

Sample Availability

Samples of the clustering structures are available from the authors.

Abbreviations

The following abbreviations are used in this manuscript:
GtACR1  Guillardia theta anion channelrhodopsin 1
ChRsChannelrhodopsins
ACRAnion channel rhodopsin
TMTransmemebranes
EVExtracellular vestibule
CCSCentral constriction site
RSBH + Protonated retinal schiff base
FTIRFourier transform infrared spectroscopy
RMSDRoot mean square deviation
ChR2Channelrhodopsin-2

References

  1. Boyden, E.S.; Zhang, F.; Bamberg, E.; Nagel, G.; Deisseroth, K. Millisecond-timescale, genetically targeted optical control of neural activity. Nat. Neurosci. 2005, 8, 1263–1268. [Google Scholar] [CrossRef]
  2. Deisseroth, K.; Feng, G.; Majewska, A.K.; Miesenböck, G.; Ting, A.; Schnitzer, M.J. Next-generation optical technologies for illuminating genetically targeted brain circuits. J. Neurosci. 2006, 26, 10380–10386. [Google Scholar] [CrossRef]
  3. Sineshchekov, O.A.; Jung, K.H.; Spudich, J.L. Two rhodopsins mediate phototaxis to low-and high-intensity light in Chlamydomonas reinhardtii. Proc. Natl. Acad. Sci. USA 2002, 99, 8689–8694. [Google Scholar] [CrossRef]
  4. Nagel, G.; Ollig, D.; Fuhrmann, M.; Kateriya, S.; Musti, A.M.; Bamberg, E.; Hegemann, P. Channelrhodopsin-1: A light-gated proton channel in green algae. Science 2002, 296, 2395–2398. [Google Scholar] [CrossRef]
  5. Nagel, G.; Szellas, T.; Huhn, W.; Kateriya, S.; Adeishvili, N.; Berthold, P.; Ollig, D.; Hegemann, P.; Bamberg, E. Channelrhodopsin-2, a directly light-gated cation-selective membrane channel. Proc. Natl. Acad. Sci. USA 2003, 100, 13940–13945. [Google Scholar] [CrossRef]
  6. Rost, B.R.; Schneider-Warme, F.; Schmitz, D.; Hegemann, P. Optogenetic tools for subcellular applications in neuroscience. Neuron 2017, 96, 572–603. [Google Scholar] [CrossRef] [PubMed]
  7. Govorunova, E.G.; Sineshchekov, O.A.; Janz, R.; Liu, X.; Spudich, J.L. Natural light-gated anion channels: A family of microbial rhodopsins for advanced optogenetics. Science 2015, 349, 647–650. [Google Scholar] [CrossRef]
  8. Sineshchekov, O.A.; Govorunova, E.G.; Li, H.; Spudich, J.L. Gating mechanisms of a natural anion channelrhodopsin. Proc. Natl. Acad. Sci. USA 2015, 112, 14236–14241. [Google Scholar] [CrossRef]
  9. Vogt, N. Algae are the best engineers of optogenetic inhibitors. Nat. Methods 2015, 12, 806–807. [Google Scholar] [CrossRef]
  10. Berndt, A.; Lee, S.Y.; Wietek, J.; Ramakrishnan, C.; Steinberg, E.E.; Rashid, A.J.; Kim, H.; Park, S.; Santoro, A.; Frankland, P.W.; et al. Structural foundations of optogenetics: Determinants of channelrhodopsin ion selectivity. Proc. Natl. Acad. Sci. USA 2016, 113, 822–829. [Google Scholar] [CrossRef] [Green Version]
  11. Li, H.; Sineshchekov, O.A.; Wu, G.; Spudich, J.L. In vitro activity of a purified natural anion channelrhodopsin. J. Biol. Chem. 2016, 291, 25319–25325. [Google Scholar] [CrossRef]
  12. Peralvárez-Marín, A.; Garriga, P. Optogenetics Comes of Age: Novel Inhibitory Light-Gated Anionic Channels Allow Efficient Silencing of Neural Function. ChemBioChem 2016, 17, 204–206. [Google Scholar] [CrossRef] [PubMed]
  13. Sineshchekov, O.A.; Li, H.; Govorunova, E.G.; Spudich, J.L. Photochemical reaction cycle transitions during anion channelrhodopsin gating. Proc. Natl. Acad. Sci. USA 2016, 113, E1993–E2000. [Google Scholar] [CrossRef] [PubMed]
  14. Yi, A.; Mamaeva, N.; Li, H.; Spudich, J.L.; Rothschild, K.J. Resonance Raman study of an anion channelrhodopsin: Effects of mutations near the retinylidene Schiff base. Biochemistry 2016, 55, 2371–2380. [Google Scholar] [CrossRef]
  15. Doi, S.; Tsukamoto, T.; Yoshizawa, S.; Sudo, Y. An inhibitory role of Arg-84 in anion channelrhodopsin-2 expressed in Escherichia coli . Sci. Rep. 2017, 7, 41879. [Google Scholar] [CrossRef] [PubMed]
  16. Yi, A.; Li, H.; Mamaeva, N.; Fernandez De Cordoba, R.E.; Lugtenburg, J.; DeGrip, W.J.; Spudich, J.L.; Rothschild, K.J. Structural changes in an anion channelrhodopsin: Formation of the K and L intermediates at 80 K. Biochemistry 2017, 56, 2197–2208. [Google Scholar] [CrossRef]
  17. Sineshchekov, O.A.; Govorunova, E.G.; Li, H.; Wang, X.; Spudich, J.L. Opposite charge movements within the photoactive site modulate two-step channel closing in GtACR1. Biophys. J. 2019, 117, 2034–2040. [Google Scholar] [CrossRef] [PubMed]
  18. Kojima, K.; Miyoshi, N.; Shibukawa, A.; Chowdhury, S.; Tsujimura, M.; Noji, T.; Ishikita, H.; Yamanaka, A.; Sudo, Y. Green-sensitive, long-lived, step-functional anion channelrhodopsin-2 variant as a high-potential neural silencing tool. J. Phys. Chem. Lett. 2020, 11, 6214–6218. [Google Scholar] [CrossRef]
  19. Kim, Y.S.; Kato, H.E.; Yamashita, K.; Ito, S.; Inoue, K.; Ramakrishnan, C.; Fenno, L.E.; Evans, K.E.; Paggi, J.M.; Dror, R.O.; et al. Crystal structure of the natural anion-conducting channelrhodopsin GtACR1. Nature 2018, 561, 343–348. [Google Scholar] [CrossRef]
  20. Li, H.; Huang, C.Y.; Govorunova, E.G.; Schafer, C.T.; Sineshchekov, O.A.; Wang, M.; Zheng, L.; Spudich, J.L. Crystal structure of a natural light-gated anion channelrhodopsin. Elife 2019, 8, e41741. [Google Scholar] [CrossRef]
  21. Dreier, M.A.; Althoff, P.; Norahan, M.J.; Tennigkeit, S.A.; El-Mashtoly, S.F.; Lübben, M.; Kötting, C.; Rudack, T.; Gerwert, K. Time-resolved spectroscopic and electrophysiological data reveal insights in the gating mechanism of anion channelrhodopsin. Commun. Biol. 2021, 4, 578. [Google Scholar] [CrossRef] [PubMed]
  22. Tsujimura, M.; Noji, T.; Saito, K.; Kojima, K.; Sudo, Y.; Ishikita, H. Mechanism of absorption wavelength shifts in anion channelrhodopsin-1 mutants. Biochim. Biophys. Acta (BBA) Bioenerg. 2021, 1862, 148349. [Google Scholar] [CrossRef]
  23. Tsujimura, M.; Kojima, K.; Kawanishi, S.; Sudo, Y.; Ishikita, H. Proton transfer pathway in anion channelrhodopsin-1. Elife 2021, 10, e72264. [Google Scholar] [CrossRef]
  24. Kuhne, J.; Eisenhauer, K.; Ritter, E.; Hegemann, P.; Gerwert, K.; Bartl, F. Early Formation of the Ion-Conducting Pore in Channelrhodopsin-2. Angew. Chem. Int. Ed. 2015, 54, 4953–4957. [Google Scholar] [CrossRef] [PubMed]
  25. Jo, S.; Kim, T.; Iyer, V.G.; Im, W. CHARMM-GUI: A web-based graphical user interface for CHARMM. J. Comput. Chem. 2008, 29, 1859–1865. [Google Scholar] [CrossRef]
  26. Lee, J.; Cheng, X.; Swails, J.M.; Yeom, M.S.; Eastman, P.K.; Lemkul, J.A.; Wei, S.; Buckner, J.; Jeong, J.C.; Qi, Y.; et al. CHARMM-GUI input generator for NAMD, GROMACS, AMBER, OpenMM, and CHARMM/OpenMM simulations using the CHARMM36 additive force field. J. Chem. Theory Comput. 2016, 12, 405–413. [Google Scholar] [CrossRef] [PubMed]
  27. Wu, E.L.; Cheng, X.; Jo, S.; Rui, H.; Song, K.C.; Dávila-Contreras, E.M.; Qi, Y.; Lee, J.; Monje-Galvan, V.; Venable, R.M.; et al. CHARMM-GUI membrane builder toward realistic biological membrane simulations. J. Comput. Chem. 2014, 35, 1997–2004. [Google Scholar] [CrossRef]
  28. Jorgensen, W.L.; Chandrasekhar, J.; Madura, J.D.; Impey, R.W.; Klein, M.L. Comparison of simple potential functions for simulating liquid water. J. Chem. Phys. 1983, 79, 926–935. [Google Scholar] [CrossRef]
  29. Phillips, J.C.; Hardy, D.J.; Maia, J.D.; Stone, J.E.; Ribeiro, J.V.; Bernardi, R.C.; Buch, R.; Fiorin, G.; Hénin, J.; Jiang, W.; et al. Scalable molecular dynamics on CPU and GPU architectures with NAMD. J. Chem. Phys. 2020, 153, 044130. [Google Scholar] [CrossRef]
  30. Darden, T.; York, D.; Pedersen, L. Particle mesh Ewald: An N· log (N) method for Ewald sums in large systems. J. Chem. Phys. 1993, 98, 10089–10092. [Google Scholar] [CrossRef] [Green Version]
  31. Hammonds, K.D.; Ryckaert, J.P. On the convergence of the SHAKE algorithm. Comput. Phys. Commun. 1991, 62, 336–351. [Google Scholar] [CrossRef]
  32. Dolinsky, T.J.; Nielsen, J.E.; McCammon, J.A.; Baker, N.A. PDB2PQR: An automated pipeline for the setup of Poisson–Boltzmann electrostatics calculations. Nucleic Acids Res. 2004, 32, W665–W667. [Google Scholar] [CrossRef] [PubMed]
  33. Isralewitz, B.; Gao, M.; Schulten, K. Steered molecular dynamics and mechanical functions of proteins. Curr. Opin. Struct. Biol. 2001, 11, 224–230. [Google Scholar] [CrossRef]
  34. Jalily Hasani, H.; Ganesan, A.; Ahmed, M.; Barakat, K.H. Effects of protein–protein interactions and ligand binding on the ion permeation in KCNQ1 potassium channel. PLoS ONE 2018, 13, e0191905. [Google Scholar] [CrossRef]
  35. Mücksch, C.; Urbassek, H.M. Accelerating steered molecular dynamics: Toward smaller velocities in forced unfolding simulations. J. Chem. Theory Comput. 2016, 12, 1380–1384. [Google Scholar] [CrossRef] [PubMed]
  36. Patel, J.S.; Berteotti, A.; Ronsisvalle, S.; Rocchia, W.; Cavalli, A. Steered molecular dynamics simulations for studying protein–ligand interaction in cyclin-dependent kinase 5. J. Chem. Inf. Model. 2014, 54, 470–480. [Google Scholar] [CrossRef]
  37. Humphrey, W.; Dalke, A.; Schulten, K. VMD: Visual molecular dynamics. J. Mol. Graph. 1996, 14, 33–38. [Google Scholar] [CrossRef]
  38. Kästner, J. Umbrella sampling. Wiley Interdiscip. Rev. Comput. Mol. Sci. 2011, 1, 932–942. [Google Scholar] [CrossRef]
  39. Akhshi, P.; Wu, G. Umbrella sampling molecular dynamics simulations reveal concerted ion movement through G-quadruplex DNA channels. Phys. Chem. Chem. Phys. 2017, 19, 11017–11025. [Google Scholar] [CrossRef]
  40. Kumar, S.; Rosenberg, J.M.; Bouzida, D.; Swendsen, R.H.; Kollman, P.A. The weighted histogram analysis method for free-energy calculations on biomolecules. I. The method. J. Comput. Chem. 1992, 13, 1011–1021. [Google Scholar] [CrossRef]
  41. Volkov, O.; Kovalev, K.; Polovinkin, V.; Borshchevskiy, V.; Bamann, C.; Astashkin, R.; Marin, E.; Popov, A.; Balandin, T.; Willbold, D.; et al. Structural insights into ion conduction by channelrhodopsin 2. Science 2017, 358, eaan8862. [Google Scholar] [CrossRef] [PubMed]
  42. Kato, H.E.; Zhang, F.; Yizhar, O.; Ramakrishnan, C.; Nishizawa, T.; Hirata, K.; Ito, J.; Aita, Y.; Tsukazaki, T.; Hayashi, S.; et al. Crystal structure of the channelrhodopsin light-gated cation channel. Nature 2012, 482, 369–374. [Google Scholar] [CrossRef]
  43. Krause, B.S.; Kaufmann, J.C.; Kuhne, J.; Vierock, J.; Huber, T.; Sakmar, T.P.; Gerwert, K.; Bartl, F.J.; Hegemann, P. Tracking pore hydration in channelrhodopsin by Site-Directed Infrared-Active azido probes. Biochemistry 2019, 58, 1275–1286. [Google Scholar] [CrossRef] [PubMed]
  44. Lórenz-Fonfría, V.A.; Bamann, C.; Resler, T.; Schlesinger, R.; Bamberg, E.; Heberle, J. Temporal evolution of helix hydration in a light-gated ion channel correlates with ion conductance. Proc. Natl. Acad. Sci. USA 2015, 112, E5796–E5804. [Google Scholar] [CrossRef] [PubMed] [Green Version]
Figure 1. The GtACR1 photocycle. The potential ion channel is indicated by the black arrow.
Figure 1. The GtACR1 photocycle. The potential ion channel is indicated by the black arrow.
Processes 11 00510 g001
Figure 2. The local conformations near CCS1 of (A) GtACR1-trans and (B) GtACR1-cis. The profiles of the potential anion conductance pathway are represented by dots through CCS1 (E68, S43, Q46, N239). (C) Comparison of distances between residue pairs of CCS1 of trans and cis conformation of GtACR1.
Figure 2. The local conformations near CCS1 of (A) GtACR1-trans and (B) GtACR1-cis. The profiles of the potential anion conductance pathway are represented by dots through CCS1 (E68, S43, Q46, N239). (C) Comparison of distances between residue pairs of CCS1 of trans and cis conformation of GtACR1.
Processes 11 00510 g002
Figure 3. Water distributions in (A) GtACR1-trans and (B) GtACR1-cis. The last 50 ns equilibrium trajectories were picked for grid water density calculation.
Figure 3. Water distributions in (A) GtACR1-trans and (B) GtACR1-cis. The last 50 ns equilibrium trajectories were picked for grid water density calculation.
Processes 11 00510 g003
Figure 4. Traction path of chloride ions in SMD simulations. The path that chloride ions are stretched at a constant speed in (A) E68p-D234p(trans) and in (B) E68p-D234p(cis). a–f means the positions that have maximum/minimum force when the Cl was pulled.
Figure 4. Traction path of chloride ions in SMD simulations. The path that chloride ions are stretched at a constant speed in (A) E68p-D234p(trans) and in (B) E68p-D234p(cis). a–f means the positions that have maximum/minimum force when the Cl was pulled.
Processes 11 00510 g004
Figure 5. PMF constructed by umbrella sampling for Cl permeation across ion channels in GtACR1-cis.
Figure 5. PMF constructed by umbrella sampling for Cl permeation across ion channels in GtACR1-cis.
Processes 11 00510 g005
Figure 6. Water density distribution in different protonation states of E68 and D234. (A) E68p-D234p, (B) E68dp-D234p, and (C) E68dp-D234dp. (D) The number of hydrogen bonds and (E) the number of water molecules in the channel. The last 50 ns equilibrium trajectories were picked for hydrogen bond and water molecule statistics.
Figure 6. Water density distribution in different protonation states of E68 and D234. (A) E68p-D234p, (B) E68dp-D234p, and (C) E68dp-D234dp. (D) The number of hydrogen bonds and (E) the number of water molecules in the channel. The last 50 ns equilibrium trajectories were picked for hydrogen bond and water molecule statistics.
Processes 11 00510 g006
Figure 7. Clustering analysis for (A) E68p-D234p, (B) E68dp-D234p, and (C) E68dp-D234dp of GtACR1. The first, second, and third columns are the representative conformations of the top three clusters in each state. The cluster percentage is displayed in each image.
Figure 7. Clustering analysis for (A) E68p-D234p, (B) E68dp-D234p, and (C) E68dp-D234dp of GtACR1. The first, second, and third columns are the representative conformations of the top three clusters in each state. The cluster percentage is displayed in each image.
Processes 11 00510 g007
Figure 8. Revised channel opening model of GtACR1.
Figure 8. Revised channel opening model of GtACR1.
Processes 11 00510 g008
Disclaimer/Publisher’s Note: The statements, opinions and data contained in all publications are solely those of the individual author(s) and contributor(s) and not of MDPI and/or the editor(s). MDPI and/or the editor(s) disclaim responsibility for any injury to people or property resulting from any ideas, methods, instructions or products referred to in the content.

Share and Cite

MDPI and ACS Style

Liu, C.; Xin, Q.; Qin, C.; Jiang, M.; Lo, G.V.; Dou, Y.; Yuan, S. The Mechanism of Channel Opening of Anion Channelrhodopsin GtACR1: A Molecular Dynamics Simulation. Processes 2023, 11, 510. https://doi.org/10.3390/pr11020510

AMA Style

Liu C, Xin Q, Qin C, Jiang M, Lo GV, Dou Y, Yuan S. The Mechanism of Channel Opening of Anion Channelrhodopsin GtACR1: A Molecular Dynamics Simulation. Processes. 2023; 11(2):510. https://doi.org/10.3390/pr11020510

Chicago/Turabian Style

Liu, Chunyan, Qi Xin, Cai Qin, Maorui Jiang, Glenn V. Lo, Yusheng Dou, and Shuai Yuan. 2023. "The Mechanism of Channel Opening of Anion Channelrhodopsin GtACR1: A Molecular Dynamics Simulation" Processes 11, no. 2: 510. https://doi.org/10.3390/pr11020510

APA Style

Liu, C., Xin, Q., Qin, C., Jiang, M., Lo, G. V., Dou, Y., & Yuan, S. (2023). The Mechanism of Channel Opening of Anion Channelrhodopsin GtACR1: A Molecular Dynamics Simulation. Processes, 11(2), 510. https://doi.org/10.3390/pr11020510

Note that from the first issue of 2016, this journal uses article numbers instead of page numbers. See further details here.

Article Metrics

Back to TopTop