On the mechanism of GIRK2 channel gating by phosphatidylinositol bisphosphate, sodium, and the Gβγ dimer

G protein–gated inwardly rectifying K+ (GIRK) channels belong to the inward-rectifier K+ (Kir) family, are abundantly expressed in the heart and the brain, and require that phosphatidylinositol bisphosphate is present so that intracellular channel-gating regulators such as Gβγ and Na+ ions can maintain the channel-open state. However, despite high-resolution structures (GIRK2) and a large number of functional studies, we do not have a coherent picture of how Gβγ and Na+ ions control gating of GIRK2 channels. Here, we utilized computational modeling and all-atom microsecond-scale molecular dynamics simulations to determine which gates are controlled by Na+ and Gβγ and how each regulator uses the channel domain movements to control gate transitions. We found that Na+ ions control the cytosolic gate of the channel through an anti-clockwise rotation, whereas Gβγ stabilizes the transmembrane gate in the open state through a rocking movement of the cytosolic domain. Both effects alter the way in which the channel interacts with phosphatidylinositol bisphosphate and thereby stabilizes the open states of the respective gates. These studies of GIRK channel dynamics present for the first time a comprehensive structural model that is consistent with the great body of literature on GIRK channel function.

ers and/or heterotetramers with each other. As their name implies, physiologically they are directly gated by G proteins, specifically the G␤␥ dimer of GTP-binding (G) proteins (G␤␥) (3)(4)(5), ones that associate with pertussis toxin-sensitive G␣ subunits. GIRK channels are dependent on phosphoinositides, especially PIP 2 , which via direct interactions stabilize their gates in the open state (6 -8). Another physiological regulator of GIRK activity is intracellular Na ϩ that is in part coordinated by a specific Asp residue found in the GIRK2 and GIRK4 channel members (9 -12). Both Na ϩ and G␤␥ have been shown to work by increasing the affinity of the channel to PIP 2 and to synergize via this mechanism in activating the channel, whereas PIP 2 itself is unable to gate efficiently the channel on its own, unlike other Kir channel family members (2). GIRK activation inhibits excitability, slowing the rate of pacemaker and atrial cell firing in the heart, while inhibiting transmitter release by presynaptic neurons or opposing excitation of postsynaptic neurons. GIRK inhibition through hydrolysis or dephosphorylation of PIP 2 causes channel desensitization that lasts as long as it takes to resynthesize PIP 2 (13,14). GIRK channels are thought to be good therapeutic targets for multiple conditions that call for a decrease in excitability or arrhythmogenesis, such as epilepsy in the brain and atrial fibrillation in the heart (15)(16)(17).
Structural studies using crystallography or computational modeling have produced 3D models of GIRK channels in complex with each and all of these physiological regulators, PIP 2 , G␤␥, and Na ϩ (4,5,8,18). Atomic resolution crystal structures of GIRK channels have been limited to a truncated GIRK2 channel construct that reproduces functional characteristics of the full-length channel (5,8,19). The GIRK2 structures have provided direct evidence for the multiple gate hypothesis of these channels. These structures point to three constrictions along the permeation pathway: the selectivity filter, where K ϩ ions are selected over other ions; the helix bundle crossing (HBC) gate, also referred in the literature as the inner helix gate, located near the inner leaflet of the plasma membrane bilayer; and a cytosolic gate that is unique to Kir channels, referred to as the G-loop gate (Fig. S1). The HBC gate is comprised by the side chains of Phe-192, whereas the G-loop gate is comprised by two methionine pairs at the top and bottom of the G loop (Fig. S2B).
The GIRK2 structures in complex with the various regulators, obtained mainly from the MacKinnon lab, have provided a wealth of structural information. However, no satisfying structural gating model has resulted that is consistent with the physiological regulation of the channel. None of the structures with physiological regulators have captured the channel gates in the open state. Even the channel bound to all its gating molecules, PIP 2 , Na ϩ , and G␤␥, which in electrophysiological experiments has been found to activate the channel the most (20), did not capture the channel gates in the open state. The structure of a low activity point mutant, GIRK2 (R201A), has come the closest to capturing the gates open, where two of the subunits bound to PIP 2 depicted open conformations, but the remaining two that were not PIP 2 bound had their gates closed (5). It was hypothesized that the crystal packing prevented PIP 2 from binding to all four subunits, leaving two of the four gates in the closed state. Based on this structure, models were created in an effort to understand what structural changes could be predicted, if all subunits would reach the open state, as expected. The binding of G␤␥ was predicted to cause an anti-clockwise movement of ϳ4°of both the cytoplasmic and transmembrane (TM) domains that put the channel in a "pre-open" state. The R201A mutant channel showed an additional ϳ4°twisting movement to cause full gate opening (5,21). Table S1 and Fig. S3 show all the crystal structures determined under different conditions with both minimal and ␣-carbon (C␣) distances between the gates of opposite-facing subunits. This heroic work has answered many questions, but at the same time has left many other questions unanswered. It has also generated puzzling conclusions, perhaps the most troubling of which is the question: does the low activity R201A mutant of the GIRK2 channel represent a valid model of a constitutively active channel that mimics the natural open channel state? Critical mechanistic questions remain, such as which gates are controlled by the two different activators (Na ϩ versus G␤␥ or both) and how the movements caused by these regulators tie into the regulation of the gates by PIP 2 . Given the nature of these dynamic questions, we pursued them using computational modeling and MD simulations based on the pre-open crystal structure of GIRK2.

Na ؉ or G␤␥ couple to distinct gates but when present simultaneously show synergism in gating the channel
Five modeled systems ranging from 160,000 to 400,000 atoms were simulated in the absence and presence of an electric field (EF) for a total of 1.0 -1.5 s, respectively (Table S2). The largest system with all the components is illustrated in Fig. S2A. K ϩ ion permeation is taken as the criterion to designate a channel conformational state under specific conditions as "closed" or "open." In the GIRK2-Apo system, no ion permeation was seen regardless of whether an EF was applied or not. In contrast, in all other systems of the GIRK2 channel, such as with PIP 2 (GIRK2-PIP 2 or GIRK2), with PIP 2 ϩ Na ϩ (GIRK2-Na ϩ ), with PIP 2 ϩ G␤␥ (GIRK2-G␤␥), and with PIP 2 ϩ G␤␥ ϩ Na ϩ (GIRK2-G␤␥-Na ϩ ), K ϩ ions permeated, albeit inefficiently, even in the presence of EF across the membrane (0.06 V/nm or V m ϭ Ϫ200 mV). However, in the GIRK2-G␤␥-Na ϩ system, an almost 10-fold greater number of ions permeated under the transmembrane EF (Table 1). The relative permeation seen in these simulations agrees with electrophysiological results (20). Moreover, in the absence of an external electric field, the MD simulations also showed an outwardly directed flow of K ϩ (intracellular to extracellular) (data not shown).
To better understand the gating mechanism, we set the production stage in the presence of an EF, because gating proved more efficient under these conditions (Table 1). All systems were equilibrated by 150 ns of simulation as seen in the rootmean-square deviation plots (Fig. S4). The average distance required to occlude permeation was 5.69 Å and was derived from the GIRK2-Na ϩ system, which has been shown experimentally to be the least conductive (20). When both the HBC and G-loop gates exceeded the minimal distance of 5.69 Å, the channel was considered to be in the open state. As designated by vertical lines in Fig. 1, this representation of activity would correspond to the open probability (P o ) of the channel in the nanosecond time scale, even though not every time the channel was permeable did an ion in fact permeate. From this simplified analysis, it became apparent that Na ϩ or G␤␥ subunits alone could open GIRK2 and could do so in a synergistic manner when applied together (Fig. 1). In fact, the synergism between the G␤␥ subunits and Na ϩ can be seen in records of single channel activity (cell-attached recordings) of the highly related GIRK4 channel. The gating behavior of certain GIRK4 mutants that mimic the Na ϩ (i.e. I229L), G␤␥ (i.e. S176P), or Na ϩ ϩ G␤␥ (i.e. I229L ϩ S176P) conditions (Fig. 2) (12, 23) bears MD simulations were run on five systems: GIRK2-Apo, GIRK2-PIP 2 (GIRK2), GIRK2-Na ϩ , GIRK2-G␤␥, and GIRK2-G␤␥-Na ϩ . Channels with their HBC and G-loop gates larger than a minimum distance (5.69 Å) that is required for permeation were considered to be open. This minimum distance estimate was derived from the least conductive system GIRK2-Na ϩ compared with the other two systems studied that have been shown experimentally to display measurable currents (i.e. GIRK2-G␤␥ and GIRK2-G␤␥-Na ϩ ). remarkable pattern similarity to the simulated records of P o , even though the time scales separating these records are at least 6 orders of magnitude apart.

Na ؉ ions predominantly stabilize the G-loop gate in the open state
The GIRK2-PIP 2 system revealed that although the HBC gate distances were distributed around the minimal distance for permeation, the G-loop gate distances became limiting (Fig. 3), explaining the inefficient open probability seen in Fig. 1. The open probability of each of the two gates individually (HBC and G loop) is shown in Fig. 4. Once Na ϩ was included in the GIRK2-PIP 2 system, there was a clear right shift of the G-loop gate minimum distributions to the open state (Fig. 3). This effect on the G-loop gate was in sharp contrast to the HBC gate, whose minimal distance distributions were not affected by Na ϩ , remaining centered around the nonpermissive 5.69 Å distance (Fig. 3). The G-loop gate P o increased by more than 7-fold, whereas the HBC gate P o , if anything, was somewhat decreased ( Table 2 and Fig. 4). Overall, the P o of both gates was simultaneously increased by ϳ6-fold ( Fig. 1 and Table 2). As a result, the HBC became the limiting gate to ion flow in the presence of intracellular Na ϩ ions.

G␤␥ dimer acts to predominantly stabilize the HBC gate in the open state
A clear right shift of the HBC gate minimum distributions to the open direction was seen when G␤␥ was included in the GIRK2-PIP 2 system (Fig. 3). The effect of G␤␥ on the G-loop gate also showed a right shift (compared with GIRK2-PIP 2 ), but there was a significant component of distance distributions that G␤␥ was not able to affect (Figs. 3 and 4 and Table 2). Opening of the G-loop gate may become a rate-limiting step to the conductive state when channels get activated by G␤␥. Because this limitation was less than that imposed by the HBC gate in the GIRK2-Na ϩ system, a higher overall P o of the chan- Unitary activity of GIRK4 * mutants (using as control the highly active GIRK4-S143T). GIRK4* mutants mimic gating by Na ϩ (I229L), G␤␥ (S176P), and G␤␥/Na ϩ (S176P and I229L) in compressed (A) and more expanded (B) time scales. Open probability of S176P (0.072 Ϯ 0.019; n ϭ 6), I229L (0.059 Ϯ 0.014; n ϭ 5), and S176P/I229L (0.322 Ϯ 0.093; n ϭ 7) are shown. These data obtained from mutants representing the endogenous gating molecules (PIP 2 , Na ϩ , G␤␥, and Na ϩ /G␤␥) are consistent with prior reports (53).  Fig. 1 were analyzed, namely GIRK2-Apo, GIRK2-PIP 2 (GIRK2), GIRK2-Na ϩ , GIRK2-G␤␥, and GIRK2-G␤␥-Na ϩ . The red dashed line indicates the cutoff distance under which no K ϩ permeation took place. GIRK2 gating mechanism by Na ؉ , G␤␥, and PIP 2 nel was achieved ( Fig. 1 and Table 2). This observation is consistent with experimental results suggesting that in comparison, gating by G␤␥ is greater than that by intracellular Na ϩ (20).

Na ؉ and G␤␥ subunits work together to activate the channel the most
When both G␤␥ and intracellular Na ϩ ions were included in the simulations, the distance distributions of both gates were predominantly in the open state ( Fig. 3), resulting in the greatest P o ( Fig. 4 and Table 2). In the presence of an EF, this condition produced the most efficient permeation of K ϩ ions (Table  1). These data are in good agreement with experimental results (Fig. 2) (20).

Rocking, tilting, and rotation movements of GIRK2 domains lead to channel activation
As discussed earlier, the only GIRK2 crystal structure depicting both the HBC and G-loop gates in the open state has been the GIRK2(R201A) ϩ PIP 2 (3SYQ) structure (8). Electrophysiological recordings of the mutant have revealed very small ionic currents (8,24), casting doubt regarding whether this mutant provides a meaningful model of the channel open structure. Principal component analysis (PCA) was performed on our simulated trajectories to discern channel domain movements linked to gating. This analysis suggests that the HBC and G-loop gates move in opposite directions as the channel gates (Movie S1). This kind of movement is in contrast to the proposal based on the 3SYQ open structure (R201A-PIP 2 ), in which the two gates move in the same direction. In fact, when we performed PCA on the crystal structures excluding the 3SYQ structure, similar results to our simulations were obtained, supporting a movement of the two gates in opposite directions (Movie S2). When PCA was repeated including the 3SYQ (GIRK2-R201A-PIP 2 ), the two gates were shown to move in the same direction as previously suggested (5). These results as a whole have suggested that the GIRK2(R201A) ϩ PIP 2 mutant involves perhaps a distinct reaction path and mechanism leading to the conformational state depicted by the crystal structure and may not represent a good model of the channel open state achieved by natural gating molecules.

Movements controlling HBC gating
Once the channel was gated from a pre-open to the open state, rocking movements of the CTD and tilting of the TM2 were observed (Movie S1). To quantitatively depict these movements two angles were defined: (a) a dihedral angle (TM2-CTD) formed by the C␣ atoms of Ser-196 with Ile-244 and Ile-281 on one hand and with Gln-287 on the other hand and (b) a TM2 tilt angle formed by the TM2 terminal and the vertical axis of the TM domain (Fig. S5A). The larger the TM2-CTD dihedral angles the more the CTD rocks outwardly, away from the vertical axis of the TM domain. Comparison of the closed and open structures reveal an outward CTD rocking movement during opening, an outward TM2 helix tilting and local twists in the G␤␥ binding ␤L-␤M and ␤D 2 -␤E 1 loops (Fig. 5). 3D plots were then employed to compare the relationship between the dihedral angle, the tilt angle, and the HBC C␣ distance. Fig. 6A shows clearly that as the TM2-CTD dihedral angles and the TM2 tilt angles increase (move outwards), so do the HBC C␣ distances. Generally speaking, increases in both the TM2-CTD dihedral and the TM2 tilt angles translate into opening the HBC gate.
It has been proposed that a four-degree anti-clockwise rotation (of the CTD relative to the TMD and as viewed from the membrane) will allow transition of the closed (3SYA or GIRK2-PIP 2 ϩ Na ϩ ) to the preopen (4KFM or GIRK2-PIP 2 ϩ G␤␥ ϩ Na ϩ ) conformation (5). We used the geometry center of the TMD and CTD excluding the N and C termini, as mentioned under "Experimental Procedures," to calculate this angle and estimated it to be ϳ3.68° (Fig. 7A), consistent with the previously proposed angle of 4°(5). However, no correlation was found between rotation angles and HBC C␣ distances in the conductive systems (GIRK2-Na ϩ , GIRK2-G␤␥, and GIRK2-G␤␥-Na ϩ ) that showed progressive opening of the HBC gate MD simulations were run on five systems: GIRK2-Apo, GIRK2-PIP 2 (GIRK2), GIRK2-Na ϩ , GIRK2-G␤␥, and GIRK2-G␤␥-Na ϩ . Channels with their HBC (A) and G-loop (B) gates larger than a minimum distance (5.69 Å) that is required for permeation were considered to be open. The minimum distance estimate was derived as in Fig. 1 from the least conductive system GIRK2-Na ϩ compared with the other two systems studied that have been shown experimentally to display measurable currents (i.e. GIRK2-G␤␥ and GIRK2-G␤␥-Na ϩ ). Table 2 The open probability (%) of the constrictions HBC, G-loop gate, and both for the last 350 ns that is in the equilibrium stage of MD in each simulation system  7B). In other words, we were unable to correlate the HBC gate opening to an anti-clockwise rotation activation mechanism alone.

Binding of two G␤␥ subunits (unlike a single Na ؉ ion) cause the cytosolic domain to rock
Na ϩ ions bind to the CD-loop of GIRK2 (and GIRK4) (10, 12) and stabilize its interactions with the G-loop gate (18). In contrast, G␤␥ binds GIRK channels at the intersubunit interface between the ␤D 2 -␤E 1 loop of one subunit and the ␤L-␤M loop of its adjacent subunit (4,5). Because the common edge between the two planes that defined the TM-CTD dihedral angle is at the ␤D 2 -␤E 1 loop level, we aimed to assess relative contributions of G␤␥ and Na ϩ to the local movements of both loops and further to the rocking motion described above.
Changes in the local conformations of the ␤D 2 -␤E 1 and ␤L-␤M loops that comprise the G␤␥-binding site were examined first. Two critical residues to the binding and activation by G␤␥, Phe-254 in the ␤D 2 -␤E 1 loop and Leu-344 in the ␤L-␤M loop (4,5,(25)(26)(27), are close to one another in the absence of G␤␥ and farther apart as G␤␥ occupies the cleft between the two loops. Fig. S6 shows that, consistent with the 4KFM crystal structure, both G␤␥-present systems (GIRK2-G␤␥/GIRK2-G␤␥-Na ϩ ) show a larger Phe-254 -Leu-344 C␣ distance compared with the systems in the absence of G␤␥ (GIRK2/GIRK2-Na ϩ ). Addition of G␤␥ that interacts with both ␤D 2 -␤E 1 and ␤L-␤M loops causes clear changes in the local conformations ( Fig. 8 and Fig. S5, C and D). The reason for the local twists was assessed by theangles as a function of the gating mole-  GIRK2 gating mechanism by Na ؉ , G␤␥, and PIP 2 cules. Significant differences in thedistribution of the systems in the presence of G␤␥ from those in its absence take place in both loops (Fig. S7). Gln-248 was reported to form contacts with Asn-75, Ser-98, and Trp-99 on the G␤␥ subunit, and mutations at Ser-98 and Trp-99 diminish G␤ activation of GIRK (28,29). Our simulation results are consistent with these experimental findings. In fact, with the exception of Ser-250, which did not display differences in thedistribution of the systems in the presence of G␤␥, all other residues have been reported to interact with G␤␥ (5). Similarly, thedistributions are consistent with the experimental work, showing that G␤␥ rather than Na ϩ could induce strong interactions with both loops except residue Ser-250. As a result, both the ␤D 2 -␤E 1 and the ␤L-␤M loops seem to act like a crank when G␤␥ subunits are present. The stronger (GIRK2-G␤␥-Na ϩ and GIRK2-G␤␥ Ͼ GIRK2-Na ϩ ) the twist rotations are, the more the CTD will rock.
For each GIRK2 subunit, only one Na ϩ binds to the ␤C-␤D 1 loop, whereas two G␤␥ subunits bind to the ␤D 2 -␤E 1 and ␤L-␤M loops, respectively (Fig. S2A). Because of its distant binding location, Na ϩ seems unable to induce strong rotations of both loops. Collectively, two G␤␥ subunits help one CTD rock more on both sides compared with the binding of a single Na ϩ ion. These differences between the extent of rocking movements caused by G␤␥ versus Na ϩ are consistent with the greater channel activation seen by the G␤␥ subunits versus Na ϩ (20).

Movements controlling the G-loop gate
We also examined the relationship of the TM2-CTD dihedral angle, the anti-clockwise rotation angle between TM and CTD, and the G-loop gate minimal distances in a 3D plot. In Fig. 6B, as the rotation angles increase and the TM2-CTD dihedral angles decrease (as shown by the arrow), the G-loop gate distances increase. In other words, an anti-clockwise rotation and/or an inwardly CTD rocking movement enlarge the G-loop gate.

PIP 2 interactions as a result of gating by Na ؉ and/or G␤␥
PIP 2 has been shown to be a necessary cofactor for activating GIRK channels (6,7). In contrast to its ability to gate other Kir channels, PIP 2 is unable to gate efficiently GIRK channels on its own but requires the presence of Na ϩ /G␤␥ (12). Furthermore, Na ϩ and G␤␥ together display synergism in activating GIRK channels (20). Thus, we proceeded to examine the interactions between residues comprising the GIRK2-binding site for PIP 2 and phosphates 4Ј (P 4 ) and 5Ј (P 5 ) atoms on its inositol ring.

Pull on the TM2 helix gets the HBC gate opened
G␤␥ relative to Na ϩ stabilizes key interactions of the positively charged Lys-194 residue with the negatively charged phosphates P 5 and P 4 of PIP 2 (Fig. 9, A and B). Meanwhile, Lys-200 in the GIRK2-G␤␥-Na ϩ system loses its interaction with PIP 2 (especially with the P 4 phosphate), stabilizing the salt bridge interactions of the Lys-199 with both phosphates of PIP 2 (Fig. 9, A and B). The CTD rocking and its anticlockwise rotation are likely to be driving the stronger binding of Lys-199 to the P 5 /P 4 atoms and the concomitant unbinding of Lys-200. This interaction rearrangement may underlie the pull on the TM2 helix to cause its tilting movement (Fig. 5) and as a consequence causing the HBC gate to open. Unlike the induction by G␤␥, Na ϩ in the GIRK2-Na ϩ system does not exert much influence on Lys-194, Lys-199, or Lys-200, and it is thus unable to influence the HBC gate (Fig. 9B). Snapshots of GIRK2 with each of its gating particles as they influence interactions with PIP 2 are also shown (PIP 2 , Fig. 9C; PIP 2 -Na ϩ , Fig. 9D; PIP 2 -G␤␥, Fig. 9E; and PIP 2 -Na ϩ -G␤␥, Fig. 9F).
The presence of PIP 2 alone with the channel in our study resulted in ion permeation, albeit low, indicating that GIRK2-PIP 2 could be active (Table 1). Although electrophysiological experiments have shown that PIP 2 on its own is not efficient to stabilize the channel gate in the open conformation (6), low basal activity can be seen (53). The low permeation in the GIRK2-PIP 2 system simulations could also be due to the preopen conformation as the initial structure for MD, which has GIRK2 gating mechanism by Na ؉ , G␤␥, and PIP 2 overcome part of the energy barrier required for the CTD rocking movement. When we took an HBC closed conformation as our initial model, for example the last snapshot of the GIRK2-Apo system, PIP 2 alone could not activate the channel even when the simulation time scale was extended to 1 s. In other words, PIP 2 alone could not overcome the CTD rocking energy barrier in this time scale. Furthermore, the HBC P o of GIRK2-PIP 2 (44.86%) seems larger than that of GIRK2-Na ϩ (35.14%) as shown in Fig. 4A and Table 2. However, in an additional 0.5-s simulation, the corresponding P o of GIRK2-PIP 2 and GIRK2-Na ϩ showed an approximate 7% decrease and 34% increase, respectively (data not shown), which suggests that PIP 2 alone is inefficient in stabilizing the HBC gate in the open state.

G-loop gate is stabilized by the interactions with the CD loop/ ␤I strand
The G-loop gate of GIRK2 consists of seven residues (MVE-ATGM). Other than the terminal methionines with longer side chains, it is difficult for other residues to induce strong shortrange van der Waals interactions. Instead, the long-range electrostatic interactions elicited by the only charged G-loop residue Glu-315 become critical for its movement. The G-loop gate has been found in a GIRK1 chimera to be stabilized in the open state by the equivalent Glu residue through a hydrogen bond of Glu-304 with His-222 of the CD loop (18). A similar hydrogen bond also exists in GIRK2. For example, in the GIRK2-Na ϩ system, the hydrogen bond between E315(N) and H233(O) is observed with significant occupancy 10.40% (Table 3). This bond is lost in GIRK2-G␤␥/GIRK2-G␤␥-Na ϩ systems because of longer distances (Fig. 10). The other G-loop gate stabilization in the closed state of the GIRK1 chimera was proposed to be the intersubunit Glu-304 -Arg-313 salt bridge (18). However, in GIRK2 this interaction seems to favor G-loop gate opening. Glu-315 (G loop)-Arg-324 (adjacent ␤I strand) was formed in the open state of the GIRK2-Na ϩ /GIRK2-G␤␥ systems (Fig. 10). During channel activation by Na ϩ , this salt bridge causes Glu-315 to be pulled outward (Fig. 11A), which facilitates the opening movement of G-loop gate via an anticlockwise rotation (Fig. 6).

PIP 2 -mediated Arg-230 -Arg-77/Asp-81 interactions affect stabilization of the G-loop gate
Simulations with the GIRK1 chimera have revealed that the N-terminal Arg-52 (Lys-64 in GIRK2)-PIP 2 salt bridge interactions in the G-loop gate stabilize the closed G-loop gate structure (18). This interaction switches to Arg-66 (Arg-77 in GIRK2)-PIP 2 to stabilize the G-loop gate open structure (18). GIRK2 gating mechanism by Na ؉ , G␤␥, and PIP 2 Figure 9. Changes in specific residue-PIP 2 interactions leading to HBC gating. A, distances between the N atom of the proposed key Lys residues (C) and the P 4 (left column) or P 5 (right column) atoms of PIP 2 taking the crystal structure (PDB code 4KFM) as a reference in red. B, PIP 2 -Lys pairwise Gibbs free energy change with a standard deviation bar. C-F, binding patterns between the key Lys that are involved in regulating HBC by the PIP 2 identified in the crystal structure (PDB code 4KFM), GIRK2-Na ϩ , GIRK2-G␤␥, and GIRK2-G␤␥-Na ϩ . The major contributors to binding are highlighted in cyan.

GIRK2 gating mechanism by Na ؉ , G␤␥, and PIP 2
Similarly, in the GIRK2-Na ϩ system open G-loop gate simulated structure, the Arg-77-PIP 2 salt bridge was established. The carbonyl oxygen atom of residue Arg-77 was positioned to interact with Arg-230, which is two positions away from the key Asp-228 residue that coordinates Na ϩ and Asp-81 in a Arg-230 -Arg-77-Asp-81 triad interaction pattern (Fig. 11A). In contrast to the GIRK2-Na ϩ system, in the GIRK2-G␤␥/ GIRK2-G␤␥-Na ϩ systems the R77-PIP 2 interactions are almost lost, which instead introduces a larger repulsion between Arg-230 and Arg-77 and consequently limits Arg-230 to interact only with residue Asp-81 (Fig. 11B). Considering the CTD anti-clockwise rotation, the outwardly pointing orientation of Arg-230 (GIRK2-Na ϩ ) is preferred to the inwardly pointing orientation (GIRK2-G␤␥/GIRK2-G␤␥-Na ϩ ), and it is more likely to stabilize the rotated conformation (Figs. 11 and  12). As the rotation angles decrease (GIRK2-Na ϩ Ͼ GIRK2-G␤␥ Ͼ GIRK2-G␤␥-Na ϩ ), the Glu-315-His-233 and Glu-315-Arg-324 stabilization is lost step by step. It is thus inferred that the PIP 2 mediated Arg-230 -Arg-77/Asp-81 interactions probably affect the stabilization of the G-loop gate in view of the CTD anti-clockwise rotation.
As a whole, the G loop is likely to be gated as follows: the larger Lys-64 -PIP 2 distance fluctuations facilitate the formation of an intersubunit His-68 (N terminus)-Val-351 (LM loop) hydrogen bond (Fig. 11 and Table 3). Based on the intersubunit Arg-77-PIP 2 interactions, Arg-230 in the CD loop interacts with Asp-81 (GIRK2-G␤␥/GIRK2-G␤␥-Na ϩ ) or Asp-81 ϩ Arg-77 (GIRK2-Na ϩ ). Different binding patterns determine whether the G-loop gate could be stabilized by the adjacent ␤I strand or the ␤I ϩ CD loop in view of CTD anti-clockwise rotation.

Discussion
In the past decade, great advances have been made in obtaining structural snapshots of a truncated GIRK2 channel that is in contact with its physiological intracellular regulators, PIP 2 , Na ϩ , and G␤␥. However, no coherent structural dynamic model has emerged that is consistent with over three decades of functional data on the gating of GIRK channels. Here, we used all-atom microsecond-scale MD simulations to discern the large movements leading to key molecular interactions that underlie channel gating by Na ϩ and G␤␥. In the conceptual model that has emerged, Na ϩ binding to the CD loop causes it to interact intramolecularly with the G-loop gate (Glu-315-His-233) and stabilize it in the open state, with the aid of an intermolecular interaction of the G loop with the adjacent ␤I strand (Glu-315-Arg-324). G␤␥ binds to the ␤D 2 -␤E 1 and ␤L-␤M cleft and causes a large change in the TM2-CTD dihedral angle, which serves as a crank to initiate a rocking movement of the CTD and, as a consequence, the tilting of TM2 helix and the opening of the HBC gate.
These movements alter the interactions of GIRK2 with PIP 2 that serve to stabilize the gates in the appropriate state causing channel gating to take place. In the case of the G-loop gate, the pivotal role of the CD loop interactions with the G loop is regulated on one hand by Asp-228 and His-233 that are utilized for Na ϩ ion coordination (8,30) and on the other hand by Arg-230 that is influenced by a slide helix PIP 2 -interacting residue, Arg- Table 3 The

key hydrogen bonds of cytosolic domain proposed to facilitate the G-loop gating
The values are given as the number/total occupancy over simulation time.  GIRK2 gating mechanism by Na ؉ , G␤␥, and PIP 2 77. Based on the intersubunit Arg-77-PIP 2 interaction, Arg-230 in the CD loop interacts with Asp-81 (GIRK2-G␤␥/ GIRK2-G␤␥-Na ϩ ) or Asp-81 ϩ Arg-77 (GIRK2-Na ϩ ), which determines whether the G-loop gate can be stabilized by the adjacent ␤I strand or ␤I ϩ CD loop (Fig. 11). On the other hand, the Arg-77-PIP 2 interaction seems to act as the connection of the G-loop gate to the HBC. The negatively charged head of PIP 2 is located between the positively charged Arg-77 and Lys-194, experiencing competitive interactions (Fig. 12). Once the Arg-77-PIP 2 interaction takes place, the Lys-194 -PIP 2 interaction is impaired (GIRK2-Na ϩ ). Otherwise, the stronger Lys-194 -PIP 2 interaction in the GIRK2-G␤␥/GIRK2-G␤␥-Na ϩ systems pulls the bottom of TM2 to open the HBC gate (Figs. 5, 9, and 12). The key amino acid residues predicted by our simulations to be involved in the stabilization of channel-PIP 2 interactions, such as Lys-194, Arg-230, Asp-81, Glu-315, and Arg-324, are highly conserved among Kir channels (Fig. S8), further underscoring their important role proposed by our study. The importance of these residues proposed by our computational studies ought to be validated experimentally in future studies. Our observations from both the HBC and G-loop gates do account for the movement of the two gates in opposite directions, as indicated by the PCA. In the presence of only Na ϩ ions, the HBC of the GIRK2-Na ϩ system is barely affected because of the limited CTD rocking movement; the G-loop gate is enlarged by an anti-clockwise rotation with the aid of the Glu-315-His-233 and Glu-315-Arg-324 interactions. On the other hand, in the presence of only G␤␥, the HBC in the GIRK2-G␤␥ system is widened by the induced larger CTD rocking; in contrast, the G-loop gate is somewhat destabilized by the loss of the Glu-315-His-233 hydrogen bond. Lastly, in the presence of Na ϩ and G␤␥, the HBC in the GIRK2-G␤␥-Na ϩ system assumes the largest size, but both the Glu-315-His-233 and Glu-315-Arg-324 stabilization interactions of the G-loop gate are lost. In other words, a larger outward CTD rocking movement enlarges the HBC but destabilizes the G-loop gate. In fact, none of the published GIRK2 WT crystal structures existed with both gates closed, based on the minimum distance (5.69 Å) required for permeation.

Hydrogen bond GIRK2-Na ϩ GIRK2-G␤␥ GIRK2-G␤␥-Na ϩ
The synergism between Na ϩ and G␤␥ in gating GIRK channels (20) was also observed in our simulations. The P o of the HBC gate is the highest in the GIRK2-G␤␥-Na ϩ system, whereas that of the G-loop gate is similar in the GIRK2-Na ϩ and GIRK2-G␤␥ systems (Table 2). Thus, the synergism is likely to correlate with the control of the HBC gate. Asp-228 acts as a key residue in the coordination of Na ϩ . Na ϩ in the GIRK2-G␤␥-Na ϩ system causes a decrease in the distance of Asp-228 -Arg-201 compared with the GIRK2-G␤␥ system (Fig. 12, C and D). A rotation by 7.74°of the C-linker is induced, which disrupts the Lys-200 -PIP 2 interaction that was present in the GIRK2-G␤␥ system. As a result, the Lys-194 -PIP 2 interaction is enhanced, and a larger and more stabilized HBC gate results in the GIRK2-G␤␥-Na ϩ system. Using both experimental and MD simulation approaches, Lys-200 was previously identified as a key player in GIRK2 channel gating by G␤␥ or ethanol (54). Disruption of the Lys-200 -PIP 2 interaction through neutralization mutations caused channel activation and enhancement of channel-PIP 2 interactions, a result consistent with the conclusions from the present simulation studies.
A recent computational study examined conduction through the GIRK2 channel pore but also gating in the absence on the GIRK2 gating mechanism by Na ؉ , G␤␥, and PIP 2 gating molecules G␤␥ or Na ϩ (55). By comparing these two systems, the authors found that PIP 2 enlarged the HBC gate, whereas the G-loop gate was found to be narrower than the HBC gate. However, the PIP 2 -bound GIRK2 channels were conductive, albeit inefficiently. These results are consistent with the results of our study. Moreover, because our study focused on gating, including the effects of the gating particles G␤␥ and Na ϩ , it elucidated the relative efficiencies of gating under different conditions (i.e. Na ϩ , G␤␥, or Na ϩ ϩ G␤␥). Our results are in complete agreement with experimental findings. Bernsteiner et al. (55) quoted the pivotal 1998 study by Huang et al. (7) as providing evidence that PIP 2 alone could activate GIRK1/4 channels. This is in fact not the case, because Na ϩ was also present in the solutions used by Huang et al. to result in efficient GIRK channel gating. Thus, the experimental literature agrees that Na ϩ or G␤␥ alone or together enhance the ability of PIP 2 to gate GIRK channels to different levels. In fact, a single point mutant has been shown to strengthen channel-PIP 2 interactions sufficiently, such that PIP 2 could now gate alone (without Na ϩ or G␤␥) (12), mimicking gating by Na ϩ (Fig. 2).
An interesting residue in the PIP 2 binding pocket is Trp-91, which also showed a dependence of its distance from Lys-194 based on the level of opening of the HBC gate. Fig. S9 shows that Trp-91 was 9.0 Å away from Lys-194 in GIRK2-G␤␥ but reached closer to Lys-194 in GIRK2-G␤␥-Na ϩ , approximately within 7-8 Å. Even though this is still not close enough for Trp-91 and Lys-194 to engage in cation-interactions, these subtle conformational differences in the Na ϩ ϩ G␤␥ system over the G␤␥ demonstrate the sensitivity of dynamic experiments (MD simulations) over static ones (X-ray crystal structures).
In our simulations, the GIRK2 system showed a pre-open to HBC-closed process, which suggested that the HBC gate closing comes ahead of the G-loop gate. Consistent with our previous results (18), it can be concluded that the open/closed movements of the HBC and G-loop gates are sequential but not simultaneous. From closed to open, the G-loop gate opening precedes that of the HBC gate (Ref. 18,GIRK1 chimera), whereas from open to close, the HBC closing precedes the G-loop gate (present study, GIRK2). Differences between GIRK2 and other GIRK homomeric (e.g. GIRK4) or heteromeric channels (e.g. GIRK1/2) remain to be explored to offer a molecular understanding of how intracellular channel modulators may differentially regulate the activity of GIRK channels.
In conclusion, the GIRK2 channel is activated by the intracellular regulators PIP 2 , Na ϩ ions, and/or G␤␥ subunits via rocking and anti-clockwise rotation movements of the CTD. Because of the lack of direct or allosteric interactions in either the ␤D 2 -␤E 1 or the ␤L-␤M loops, the Na ϩ ion action is limited to the control of the G-loop gate through the same subunit Glu-315-His-233 and adjacent subunit Glu-315-Arg-324 hydrogen bonds, thus inducing an anti-clockwise rotation. In contrast, binding of the G␤␥ subunits in the cleft created by the ␤D 2 -␤E 1 and the ␤L-␤M loops of adjacent subunits twists both the loops in a crank-like manner to rock the CTD and pull the HBC gate to open with the aid of Lys-194 -PIP 2 salt bridge interactions. The interactions of Arg-77 (slide helix)-PIP 2 seems to be what connects the G-loop to the HBC gate. Na ϩ enlarges the G-loop gate and drives formation of the Arg-77-PIP 2 interaction that in turn destabilizes the Lys-194 -PIP 2 interaction and the conductive state of HBC. The synergism of Na ϩ ϩ G␤␥ is likely to be attributed to a 7.74°rotation of the C-linker that as a result disrupts the Lys-200 -PIP 2 and enhances the Lys-194 -PIP 2 attraction.

Setting up GIRK2 model systems for MD simulations
The crystal structure of GIKR2 channel (PDB code 4KFM), which is considered to be in the pre-open state, was used to build the initial models by adding different endogenous regulators, such as PIP 2 , Na ϩ , and G␤␥. To characterize the role of each regulator alone or in combination, five simulation systems in total were constructed, which are termed GIRK2, GIRK2-PIP 2 , GIRK2-PIP 2 -Na ϩ (GIRK2-Na ϩ ), GIRK2-PIP 2 -G␤␥ (GIRK2-G␤␥), and GIRK2-PIP 2 -G␤␥-Na ϩ (GIRK2-G␤␥-Na ϩ ), respectively (Table S2). All the model systems were built by taking the following steps: (a) All the missing side chains of residues, Ile-55, Arg-73, Glu-127, Phe-141, Lys-165, Lys-301, and Glu-303 in GIRK2; Arg-42 and Arg-214 in the G␤ subunit; and Glu-58 and Glu-63-Phe-67 in the G␥ subunit (PDB code 4KFM) were added by Discovery Studio 2017 software. (b) All the hydrogen atoms were added using the Hϩϩ website server (http://biophysics.cs.vt.edu/) 5 (31,32). The protonated states of the titratable residues were determined by pK a calculations at neutral physiological conditions (pH 7.0). (c) The PIP 2 molecule was built based on the DiC1PIP 2 in the crystal structure (PDB code 4KFM) by adding alkyl tails using Maestro (Schrödinger, LLC). The geometry of the PIP 2 structure was optimized by ab initio quantum chemistry at the Hartree-Fock/6 -31G* level, followed by restrained electrostatic potential charge calculations using the GAUSSIAN 09 program (33). The antechamber module of AmberTools17 was employed to generate the required force field parameters for PIP 2 based on the derived restrained electrostatic potential charges and the GAFF2 force field (34 -36). (d) The GIRK2 channel and complexes were immersed in explicit lipid bilayer of 1-palmitoyl-2-oleoyl-sn-glycero-3-phosphocholine, 1-palmitoyl-2oleoyl-sn-glycero-3-phosphatidylethanolamine, 1-palmitoyl-2oleoyl-sn-glycero-3-phosphatidylserine, and cholesterol with molecular ratio of 25:5:5:1 (37) and a water box (116 ϫ 116 ϫ 150 Å 3 for GIRK2, GIRK2-PIP 2 , and GIRK2-Na ϩ and 162 ϫ 162 ϫ 170 Å 3 for GIRK2-G␤␥ and GIRK2-G␤␥-Na ϩ , a box edge of at least 14 Å from the protein periphery in each dimension using periodic boundary conditions) by using the CHARMM-GUI Membrane Builder webserver (http://www.charmm-gui.org/ ?docϭinput/membrane) 5 (38 -41). 150 mM KCl was added into the system with the K ϩ ions, and water molecules from the crystal structure retained. (e) The tleap module of Ambertools17 was used to neutralize the complexes by adding additional K ϩ or Cl Ϫ ions prior to generation of the topology and coordinates files. The FF14SB, LIPID17, and GAFF2 force fields were chosen for protein, mixed lipid membrane, and PIP 2 , respectively. The simulation systems obtained contain a range from 160,000 to 400,000 atoms (Table S2). 5 Please note that the JBC is not responsible for the long-term archiving and maintenance of this site or any other third party hosted site.

All-atom microsecond-scale MD simulations
A two-stage energy minimization protocol, steepest descent algorithm (10,000 steps) and conjugate gradient (10,000 steps) was performed for each model system. The systems were heated from 0 to 300 K using Langevin thermostat algorithm with a 0.5-fs time step to avert internal disturbance. In the heating stage, the protein and lipid bilayer were initially fixed to remove any potential steric clashes from K ϩ or Cl Ϫ ions, and water molecules; followed by gradually reduced position restraints on the protein and membrane (10 to 0.1 kcal/mol⅐Å 2 in 10 steps of total 20 ns). 0.5-s MD simulations were conducted using the constanttemperature, constant-pressure ensemble (NPT) without electric field, followed by the same time scale using the constant-temperature, constant-volume ensemble (NVT). To accelerate the permeation events in the limited time scale of simulations, an external voltage 0.06 V/nm (18,42) was employed, under which the secondary structures of all four subunits of GIRK2 were well-maintained. However, higher electric fields resulted in structural instability. The PMEM-D.CUDA program in AMBER16 was used to conduct the simulations. Long-range electrostatics were calculated using the particle mesh Ewald method with a 10-Å cutoff. A 4-fs time step by employing hydrogen mass repartition algorithm for system solutes (43) was used to accelerate the MD simulations. The SHAKE algorithm was applied on the solvent molecules.

Analysis of MD simulation results
All the analysis was done on the trajectories with external electric field, unless otherwise mentioned. Geometry analysis (distances and dihedral angles), PCA, and generalized Born surface area binding free-energy calculations were implemented by Amber16 and Ambertools17 (44). The channel-gating mechanism was studied by PCA based on the concatenated trajectories of GIRK2 (pre-open to closed) and GIRK2-G␤␥-Na ϩ (pre-open to open). For comparison, the PCA was also performed on available GIRK2 crystal structures (PDB codes 3SYA, 3SYC, 3SYO, 3SYP, and 4KFM) with and without the "open" conformation R201A with PIP 2 (PDB code 3SYQ). To evaluate the TM-CTD rotation angle, larger-displacement residues in the N and C termini were excluded. To make the angle more reasonable, the relative tilting of TM domain to CTD occurred during simulations was removed by our program in Python 2.7.5, implemented in MDAnalysis 0.18 (45,46). The delta TM-CTD rotation angle was obtained by taking the same reference, a closed conformation (PDB code 3SYA) (5). The HOLE program version 2.2 was used for analysis of the dimensions of the pore in GIRK2 channel structures (47,48). The sequence conservation analysis was performed using the Con-Surf server (http://consurf.tau.ac.il/) 5 (49 -51). VMD 1.9.3 was used for structure rendering (52). All the calculations were conducted on servers equipped with NVIDIA Tesla K80, K40, and GeForce GTX-1080 graphical cards. The executables were built under the Community Enterprise Operating System (CentOS) version 7 and NVIDIA Compute Unified Device Architecture (CUDA) toolkit version 8.0.

Ion channel expression
Various ion channel cDNAs were generated by introducing point mutations on the background of the channel GIRK4 (S143T) (56) using the QuikChange site-directed mutagenesis kit (Stratagene, La Jolla, CA). cRNAs were in turn generated by in vitro transcription using the Message Machine kit (Ambion, Austin, TX). The cRNA concentrations were estimated by comparing the fluorescence intensity of diluted cRNA with cRNA marker (Gibco) which were electrophoresed in parallel on formaldehyde gels. Ion channel proteins were expressed in Xenopus oocytes by injecting the desired amount of cRNA into the oocytes. The desirable expression level for single channel recording was achieved by incubating the oocytes at 18°C for 1-2 days after injecting 0.01-0.1 ng cRNA/oocyte. Oocytes were isolated and microinjected as previously described (22). The institutional animal care and use committee-approved protocol for use of Xenopus laevis frogs to isolate oocytes at Northeastern University was last approved on December 2018 by protocol 17-0102R-A1.

Single-channel recording and analysis
The single-channel activity was recorded using an Axopatch 200A amplifier (Axon Instruments). The pipette solution contained 96 mM KCl, 1 mM MgCl 2 , and 10 mM HEPES (pH 7.40). The bath solution contained 96 mM KCl, 5 mM EGTA, and 10 mM HEPES (pH 7.40). 100 M gadolinium was routinely included in the pipette solution to suppress native stretch channel activity in the oocyte membrane. Chemicals were purchased from Sigma. Single-channel currents were filtered at 1-2 kHz with a six-pole low-pass Bessel filter, sampled at 5-10 kHz, and stored directly into the computer's hard disk through the DIGI-DATA 1200 interface (Axon Instruments). Single-channel analysis was carried out with pClamp8 (Axon Instruments). The open probability of a single-channel in a recording was calculated by dividing the sum of the durations that the channel dwelling at the open state by the total length of recording. The length of the recordings used for computing the channel open probability ranged from 100.0 to 1804.8 s (573.8 Ϯ 499.9 s, n ϭ 18).