Thermodynamic coupling between activation and inactivation gating in potassium channels revealed by free energy molecular dynamics simulations

K⁺ channels control the flow of K⁺ ions across cell membranes principally via two different mechanisms involving structural rearrangements either at the intracellular gate or at the selectivity filter. The conformational states of these two elements are dynamically coupled. Opening of the intracellular gate via external stimuli (activation) results in a transient period of sustained ion conduction until the selectivity filter undergoes a spontaneous conformational change (inactivation) toward a collapsed or constricted nonconductive state (Hoshi et al., 1991; Yellen, 1998; Gao et al., 2005; Cordero-Morales et al., 2006a, 2007; Cuello et al., 2010b,a). Activation and inactivation gating are allosterically coupled: an open gate favors an inactivated filter, and a closed gate stabilizes a conductive filter (Panyi and Deutsch, 2006; Cuello et al., 2010a). Ultimately, the rates of entry into and recovery from inactivation govern the total amount of K⁺ membrane permeation available in response to a given external stimulus. The complex interplay between activation and inactivation, a hallmark of all K⁺ channels, underlies a wide range of electrical events in excitable and nonexcitable cells. Although it is customary to discuss activation and inactivation of K⁺ channels in terms of two functional “gates,” the molecular basis of how those structural elements might be allosterically coupled is not well understood.

The bacterial KcsA K⁺ channel from Streptomyces lividans provides an excellent prototypical model system to understand activation–inactivation coupling. Although the KcsA channel comprises only a pore domain without the complex machinery of voltage-gated K⁺ channels, it contains all the molecular elements necessary for ion conduction, activation, and inactivation gating. In the closed state, the inner helices of each of the four subunits constituting the pore domain adopt a straight configuration, forming a tight bundle that occludes the permeation pathway at the intracellular end (Doyle et al., 1998; Zhou et al., 2001). Opening of the intracellular gate is associated with a hinge-bending motion that allows these four helices to splay outward, creating a wide aqueous vestibular entryway leading up to the narrow selectivity filter (Perozo et al, 1999; Liu et al., 2001; Jiang et al., 2002; Cuello et al., 2010b).

A recent analysis based on multiple x-ray structures of open and partially open KcsA revealed the mechanism by which movements of the inner bundle activation gate is linked to conformational changes within the selectivity filter (Cuello et al., 2010a,b). The analysis highlighted the important role of Phe103, a residue located along the inner helix near the hinging position associated with the opening of the intracellular gate. In the open state, it was observed that the aromatic side chain of Phe103 changes its rotameric state, leading to an apparent steric clash with residues of the selectivity filter from the neighboring subunit (Fig. 1). Mutations at position 103 were shown to affect gating kinetics in a size-dependent way, with small side-chain substitutions leading to considerably slower entry into inactivation. Strikingly, an alanine residue substituted at position 103 abolished inactivation in KcsA (Cuello et al., 2010a). Average energies calculated from equilibrium molecular dynamics (MD) simulations indicated that, among the different residues along the inner helix, the strongest energetic interactions with the selectivity filter involved the side chain of Phe103 (Cuello et al., 2010a). These observations led to the proposition that Phe103 is the structural link responsible for the allosteric communication between the intracellular gate formed by the inner helical bundle and the selectivity filter.

Although this analysis based on static x-ray structures and equilibrium MD trajectories provided a plausible and coherent picture, the functional and crystallographic results beg for a more quantitative characterization of the nature of the microscopic interactions involving the residue at position 103 with the rest of the protein structure, and how it underlies the allosteric mechanism coupling the activation gate to the selectivity filter. In particular, knowledge of average interaction energies as calculated by MD is insufficient to explain the factors underlying the structural stability of functional states. The latter is governed by relative free energies, and this is the critical information that is ultimately needed for a complete and quantitative analysis of the wild-type (WT) and mutant channels.

Here, alchemical free energy perturbation (FEP) MD simulations (FEP/MD) of the F103A mutation in KcsA are used to elucidate the thermodynamic basis for the coupling between the intracellular gate and the selectivity filter. The results of the calculations are in good agreement with experiments. This then enables us to carry out a detailed dissection of the microscopic interactions underlying the relative stability of the two functional states of the filter. Decomposition of the free energy differences in terms of repulsive, van der Waals, and electrostatic contributions reveals the energetic basis underlying the absence of inactivation in the F103A mutation in KcsA. When the intracellular gate is open, the large aromatic phenylalanine side chain at position 103 appears to destabilize the conductive conformation of the selectivity filter because of an unfavorable steric van der Waals clash with Ile100 in an adjacent subunit. This idea is verified experimentally by examining the behavior of the I00A KcsA mutant, thus lending support to a mechanical coupling mechanism of activation and inactivation that is highly cooperative in nature.

Experimental macroscopic currents through the KcsA channel were recorded in response to an activating pH jump (from 8 to 3) using the liposome patch-clamping methodology that was described previously (Cordero-Morales et al., 2006a). For KcsA-WT and F103A mutant channels, macroscopic currents were obtained at +100 and −100 mV from KcsA-reconstituted liposomes in 200 mM of symmetrical conditions.

Adiabatic energy maps were calculated from scanning the rotameric energy surface along the χ₁ (N-CA-CB-CG) and χ₂ (CA-CB-CG-CD1) dihedral angles of Phe103 in one of the four subunits of KcsA in its closed (Protein Data Bank accession no. 1K4C) and open (Protein Data Bank accession no. 3F5W) conformations. For the sake of simplicity, the interaction was computed by including only a single subunit of the KcsA channel in vacuum. Atoms at distances greater than 10 Å from the center of mass of Phe103 were fixed in space; all other non-hydrogen atoms, except for those belonging to Phe103, were harmonically restrained with a force constant of 100 kcal/mol/Ų. The values of χ₁ and χ₂ were systematically varied along a grid, and new conformations of the Phe103 side chain were built using the internal coordinates module of the CHARMM program (Brooks et al., 2009). At every grid point, the angles were restrained with a force constant of 1,000 kcal/mol, the structure was minimized with 1,000 steps of steepest descent and 500 steps of adopted basis Newton-Raphson minimization, and the final energy was recorded without the dihedral restraint.

The x-ray structure of the KcsA K⁺ channel in the open-inactivated (oi) state with the intracellular gate open by 32 Å, as measured from the distance between the Cαs of Thr112 of opposing subunits (Protein Data Bank accession no. 3F5W), was used (Cuello et al., 2010b). Missing atoms at the N and C termini of 3F5W were modeled from the x-ray structure of KcsA in the closed state (Protein Data Bank accession no. 1K4C) after local alignment (Zhou et al., 2001). The structure of the channel with an open intracellular gate and open-conductive filter (oo) was modeled by replacing residues 70–80 of the selectivity filter with the coordinates from the x-ray structure with a closed intracellular gate and an open-conductive filter (Protein Data Bank accession no. 1K4C) (Zhou et al., 2001). In both cases, the simulated model comprises the channel structure embedded in a dipalmitoylphosphatidylcholine (DPPC) bilayer surrounded by an aqueous KCl salt solution. The models contained the KcsA tetramer (404 amino acid residues), 112 DPPC molecules, 7,325 water molecules, and 2 K⁺ ions in the pore at positions S1 and S4, and S2 and S4, for inactive and conductive filters, respectively. To ensure electrical neutrality and mimic a 150-mM KCl concentration, 18 K⁺ and 32 Cl⁻ ions were added to the bulk solution. The glutamic acid at position 71 was protonated (Bernèche and Roux, 2002), and the termini were patched with standard N- and C-terminal capping groups. For stage 1 of the FEP calculations, an alanine residue was substituted at position 103 in all four KcsA subunits and carefully equilibrated. All simulation systems were set up using the CHARMM program (Brooks et al., 2009) according to a methodology described previously (Bernèche and Roux, 2000).

The free energy for alchemically mutating an alanine into a phenylalanine at position 103 was calculated as a step-by-step reversible work. The alchemical transformation was performed in each of the four subunits of the KcsA channel, and thus, the FEP/MD computation results in the total observable relative free energy caused by a F103A mutation inserted into a tetrameric channel. Each step consists of an independent MD simulation corresponding to a different microscopic state, which allows the computation to be run efficiently in parallel. The FEP/MD simulations were performed on a fully solvated KcsA system in a lipid bilayer with ions as described above. Spurious instabilities can occur when simulating a moderate to low-resolution structure, which might take a long time to equilibrate, thus causing long-time convergence problems in the FEP/MD calculations. In view of the moderate resolution (3.2 Å) of the x-ray structure for the open-inactivated state (Cuello et al., 2010b), a weak harmonic restraint of 0.5 kcal/mol/Ų was applied to non-hydrogen atoms of the backbone at a distance of 10 Å away from the four Phe103 residues during the FEP/MD calculations to avoid deviations away from the open-channel structure determined by x-ray crystallography (Cuello et al., 2010b). This is justified in the present case because we aim to extract the most meaningful information from the FEP/MD simulations in a reasonable amount of computational time.

The calculation involved two major stages, depicted in Fig. 2. In the first stage, the FEP/MD simulations were performed with only an alanine at position 103, without the presence of a virtual aromatic group attached to the carbon β atom (CB) and decoupled from the surroundings. The methyl group (CH3) of the alanine was transformed into a methylene group (CH2). The second stage involves materializing the remainder of the phenylalanine side chain. This part of the calculation was performed with only a phenylalanine at position 103, without the presence of a virtual HB3 atom attached to the CB and decoupled from the surroundings. In practice, these non-interacting atoms would not affect the calculated relative free energies and can be left out for the sake of simplicity. Because the same FEP/MD procedure was followed for both the oo (open inner gate with open-conductive filter) and oi (open inner gate with inactivated filter) states of the channel leading to appropriate endpoints, the procedure leads to a properly closed thermodynamic cycle. In step 1a, a harmonic dihedral restraint with a force constant of 20 kcal/mol/rad² was applied on the χ₁ of ala103. In step 1b, the atom type of CB was changed from a methyl group carbon (CT3) to that of a methylene group atom (CT2), and the Lennard-Jones (LJ) 6–12 interactions of HB3 were switched off as its charge was transferred to the CB in the presence of the χ₁ dihedral restraint.

In steps 2a, 2b, and 2c, the repulsive, van der Waals dispersive, and electrostatic contributions to the nonbonded interactions of the aromatic group with its surroundings were progressively switched on in the presence of harmonic dihedral restraint potentials acting on the χ₁ (N-CA-CB-CG) and χ₂ (CA-CB-CG-CD1) dihedral angles of the phenylalanine side chain.

The dihedral restraint was applied with a force constant of 20 kcal/mol/rad², using the x-ray structure as a reference. For steps 2a and 2b, the LJ 6–12 potential was separated into purely repulsive and dispersive components using the Weeks-Chandler-Andersen (WCA) prescription (Weeks et al., 1971; see also Deng and Roux, 2004). Specifically, the core repulsion diverges as r^-12 at short distance, whereas the van der Waals dispersion interaction goes as r^-6 at long distance. The WCA separation of the LJ 6–12 potential is a convenient method to construct these two components based on the first derivative (i.e., force), which simplifies the interpretation of the free energy decomposition (Deng and Roux, 2004). Lastly, all dihedral restraint potentials were removed in the final step, 2d. Such dihedral restraints help improve the sampling and convergence by reducing the amount of fluctuations during intermediate stages of the FEP/MD calculation for steps 2a, 2b, and 2c, and do not affect the final relative free energies.

The progression of the repulsive, dispersive, and electrostatic contributions in the oo and oi states (steps 2a, 2b, and 2c) are shown in Fig. 5. For each stage, two independent simulations were run for the initial and final states of each window, and the total simulation time was 100 ps. Data from the first 10 ps was discarded for equilibration. All the FEP/MD calculations were performed using the PERT module of the CHARMM program (Brooks et al., 2009). The final results are given in Table II, with a complete breakdown of the various contributions to the total free energy.

The recent crystallographic data on KcsA shows that the opening of the intracellular gate is coupled to a conformational change taking place within the selectivity filter toward a constricted nonconductive state (Cuello et al., 2010b). Furthermore, a change in the rotameric state of the aromatic side chain of Phe103 is also observed. The closed (1K4C) and fully open channel (3F5W) structures with the highlighted Phe103 are shown in Fig. 1. By virtue of its location, Phe103 is an obvious candidate residue that could relay information about the conformation of the intracellular gate to the selectivity filter. It is located along the inner helix, in the vicinity of both the selectivity filter and the helix hinge bending associated with the opening of the intracellular gate. It is unclear, however, why the rotameric state of Phe103 changes upon the opening of the activation gate.

To gain more insight into the microscopic factors responsible for the change in the rotameric state of Phe103, we compared the two-dimensional maps of the energy surfaces along the χ₁, χ₂ dihedral angles of the Phe103 side chain in the closed and open state of KcsA. We emphasize that, although these are only adiabatic (energy-minimized) maps and not free energy surfaces, the results are expected to be meaningful and informative because the allowed conformations of the side chain are dominated by local steric effects. Fig. 3 shows that that there is a shift in the position of the energy minimum of the side chain when the intracellular gate opens. For the purpose of this study, the open-inactivated x-ray structure of the KcsA channel with its intracellular gate opened by 32 Å (Protein Data Bank accession no. 3F5W) is used. When the intracellular gate opens, the energetically favorable rotamer of the Phe103 side chain changes from an “up” state (Fig. 3, A and B, white square: χ₁ = −68; χ₂ = 105), in which the aromatic ring is tucked along the inner helix in close proximity to the selectivity filter, to a “down” state (white triangle: χ₁ = −155; χ₂ = 41), in which the aromatic ring points away from the helix. The change reflects the preferred rotamers commonly seen in the crystal structures of the closed (1K4C) and open (3F5W) states (see Fig. 1 A). An interesting exception, however, is the synthetic KcsA with a d-alanine substituted at the position of Gly77 (Protein Data Bank accession no. 2IH1). This is the single occurrence of a crystal structure corresponding to a KcsA with its intracellular gate closed and an open-conductive selectivity filter with the ring of Phe103 pointing down (Valiyaveetil et al., 2006). This suggests that both rotamers of Phe103 are allowed when the intracellular gate is closed, although the up state appears to be more frequently observed. This agrees with the computed adiabatic map of KcsA with its intracellular gate closed (Fig. 3 A); although both the up and down states of F103 are energetically allowed, the up state is more favorable.

The energy maps confirm that these two rotameric states are favored in their respective conformations and that the down state is highly favored when the intracellular gate is open. The change in rotameric state is essentially driven by the conformation of the backbone associated with channel gating. In the two-dimensional maps, the up rotameric state observed in the closed channel (Fig. 3 A, white square) corresponds to an unfavorable energy region in the open channel (Fig. 3 B, white square). As shown in Table I, the distance between the bulky residues of the inner helix (Ile100 and Phe103) and those located near the tip of the selectivity filter (e.g., Ala73 and Thr74) decreases significantly, by ∼1.5 to 2.0 Å, in going from the closed state to the open-inactivated state. Paradoxically, although the inner helices bend and splay apart upon the opening of the intracellular activation gate, there is actually a slight compression of the upper region of the inner helix (above the sharp bend) against the bottom of the selectivity filter. This is shown in Table I; there is a reduction of the distance between the bulky residues Ile100 and Phe103 of the inner helix and Ala73, Thr74, and Thr75 at the bottom of the selectivity filter. These backbone displacements have a strong impact on the accessible rotamers of Phe103. Examination of the conformations corresponding to the high energy region in Fig. 3 shows that the side chain of Phe103 cannot remain in the same rotameric state (χ₁ = −68; χ₂ = 105) because of clashes with Ala73 and Thr74 (junction between the pore helix and selectivity filter), and Leu36 (outer helix) (Fig. 3 C). In the closed state (1K4C), the Cα-Cα distance between Phe103 and Ala73 is 9.2 Å, and it is 8.3 Å between Phe103 and Thr74. In contrast, those distances are reduced to 6.8 and 6.9 Å in the open state. As a result, the rotameric state of Phe103 in the closed state becomes energetically forbidden due to the steric clash between the aromatic side chain and some residues within the same subunit. Thus, in response to outward bending of the inner helix associated with the opening of the intracellular activation gate, the side chain of Phe103 is essentially obligated to change its rotameric state.

By virtue of its response to the intracellular gate and its spatial proximity to the selectivity filter, Phe103 represents an important structural link that couples channel opening and C-type inactivation. To transmit information from the intracellular gate to the filter, however, the change in the Phe103 rotamer needs to affect the stability of the conductive pore. To help formulate this issue more quantitatively, it is useful to consider the minimal number of channel states that could possibly result from these two dynamical elements. The scheme involves transitions between four states: closed inner gate with open-conductive filter (co), open inner gate with open-conductive filter (oo), open inner gate with inactivated filter (oi), and closed inner gate with inactivated filter (ci). In this way, the observed kinetic behavior of WT and mutant KcsA can be described in terms of equilibrium population shifts among the four preexisting states (co, oo, oi, and ci) (Yifrach et al., 2009; Cuello et al., 2010a). The effectiveness of the mutation to alter entry into inactivation of the KcsA channel once the intracellular gate is open is governed by the probabilities P_(oo) and P_(oi). It is advantageous to address the issue from a thermodynamic standpoint by considering the relative free energy of the oo and oi states,

Such an expression also holds for mutants of the channel.

The experimentally observed effects of substituting the phenylalanine at position 103 by an alanine observed experimentally are shown in Fig. 4. The macroscopic current response upon a pH jump shows that the WT channel inactivates fairly rapidly, whereas the mutant F103A remains predominantly conductive. The mutation F103A sharply reduces the rate and extent of C-type inactivation. Qualitatively, inactivation appears to be essentially abolished by the F103A mutation. However, it is necessary to characterize this observation more quantitatively to compare with the results of the free energy simulations. The most reliable measure of the equilibrium effect of the mutation is via the I/I_max ratio. This allows us to estimate P_(oo), the steady-state probability of the channel to have its filter in the conductive state while the inner activation gate is open. For the sake of simplicity, we consider the behavior of the channel under a positive applied potential (outward currents) in the present analysis. As observed in Fig. 4, the inactivation of the KcsA selectivity filter is partly voltage dependent. At negative potential, the channel inactivates more rapidly because there are additional processes associated with the side chain Glu71 (Cordero-Morales et al., 2006b). It is sufficient at this point to limit the analysis to the outward currents. Based on Fig. 4, the selectivity filter of KcsA-WT has a probability of ∼0.05 on average to be in the conductive state when the inner activation gate is open; P_(oo) ranges from 0.02 to 0.07 based on the analysis of three traces (not depicted). In contrast, the selectivity filter of the F103A KcsA has a probability of ∼0.90 to be in the conductive state when the inner activation gate is open; P_(oo) ranges from 0.86 to 0.90 based on the analysis of five traces (not depicted). For the WT channel, these results imply that the probability ratio of the oo state relative to the oi state is ∼0.05/0.95 at positive voltage, corresponding to a ΔG_(oo),(oi) that is roughly equal to +1.8 kcal/mol for the WT channel. For the F103A mutant, the probability ratio of the oo state relative to the oi state increases to ∼0.90/0.10, yielding a free energy ΔG_(oo),(oi) of roughly −1.3 kcal/mol. Combining the results from several macroscopic current traces for WT KcsA and the F103A mutant (not depicted), the net effect of the mutation on the relative stability of the oo and oi state, ΔΔG, is found to vary between −3.6 and −2.6 kcal/mol (the number is negative because the mutation makes the state oo more stable relative to the state oi).

Computations can be used to gain more insight into the role of the F103 side chain and the coupling between the activation gate and the selectivity filter. A simple “brute force” MD simulation approach would seek to generate long trajectories to observe directly the relative probabilities P_(oo) and P_(oi) for the WT and mutant channel. Because inactivation in KcsA is a slow process taking place over hundreds of milliseconds, however, such an approach would not be practical. Alternatively, the problem can be treated thermodynamically by considering the relative free energy of the different functional states. Accordingly, the observable effect of a mutation on the relative probabilities of the oo and oi state can be expressed as,

Quantities such as $ΔG(oo)F103A,WT=[G(oo)F103A−G(oo)WT]$ and $ΔG(oi)F103A,WT=[G(oi)F103A−G(oi)WT]$ represent free energy differences between the WT and mutant channel in a given conformational state (oo or oi).

Although experimentally inaccessible, such free energy differences can be calculated via alchemical FEP/MD simulations, in which the step-by-step reversible thermodynamic work for mutating an amino acid is calculated. Each stage consists of an independent MD simulation with a different set of options and parameters. The physically intuitive route would compare the change in the probabilities for inactivation in the WT and mutant (F103A). The alchemical FEP/MD simulation approach involves comparing the difference in free energy for mutating (in silico) the alanine to a phenylalanine in the inactive and open-conductive states.

The results of the FEP/MD simulations are given in Table II. The total free energy difference to replace the phenylalanine at position 103 by an alanine in KcsA is consistent with a stabilization of the conductive state relative to the inactive state by −5.3 kcal/mol. This result is consistent with the important observation that the F103A mutation relieves inactivation (Cuello et al., 2010a). Moreover, the calculated ΔΔG is also in good accord with the relative free energy estimated from macroscopic current experiments, which ranges from −3.6 to −2.6 kcal/mol (see above). Once it is established that the computational result for the F103A mutation is consistent with experiment, the FEP/MD calculations give us an opportunity to gain a deeper insight into the molecular interactions driving channel inactivation. This can be achieved by considering the various free energy components of the step-by-step alchemical F103A transformation used in the FEP/MD simulations.

A particularly useful framework for dissecting the interaction of a side chain with its surroundings in terms of repulsion, dispersion, and electrostatics can be formulated based on methods developed for the calculation of the solvation of molecular solutes (Deng and Roux, 2004). The decomposition proceeds from the general observation that all intermolecular forces are dominated by harsh short-range repulsive interactions arising from the Pauli exclusion principle, long-range van der Waals attraction as a result of quantum mechanical dispersion, and electrostatic interactions from a non-uniform molecular charge distribution. In standard biomolecular potential functions, electrostatic interactions are represented on the basis of Coulomb interactions between partial charges, whereas nonpolar interactions are modeled with LJ 6–12 potentials. The latter can be separated into purely repulsive and attractive parts according to the WCA prescription (Weeks et al., 1971; Deng and Roux, 2004). In the following, the repulsive and dispersive contributions are based on the WCA separation of the LJ potential. If the repulsive, dispersive, and electrostatic parts of the interactions between the side chain and its surroundings are thermodynamically switched on in a step-by-step fashion using an FEP/MD simulation, the free energy difference arising from the nonbonded interactions may be written as, ΔΔG_int = ΔΔG_rep + ΔΔG_dis + ΔΔG_elec, where

and

Such a free energy decomposition must be interpreted with caution, however, as the various terms are conditionally defined by the order in which the interactions are introduced in the FEP/MD simulations (see Materials and methods). Fundamentally, any thermodynamic decomposition of a free energy in a strongly coupled system is a theoretical construct that is path dependent. The step-by-step FEP/MD method that we use here consists in turning on the repulsive interaction, then the van der Waals dispersion, and then the electrostatics. This path makes it possible to attribute the free energy contributions to a specific interaction, but it is dependent on the sequence chosen. For example, the dispersion is calculated once the repulsion is already active, and the electrostatics is calculated once both the repulsion and the dispersion are already active. Changing the order of how these interactions are turned on, however, could lead to nonsensical results. For example, it is not possible to turn on the dispersion or the electrostatics before the core repulsion is already active because the energy becomes divergent. More importantly, the step-by-step FEP/MD decomposition method can provide useful insight, as demonstrated by previous applications to the solvation free energy of simple compounds (Deng and Roux, 2004), and on the binding of ligands to proteins (Deng and Roux, 2009).

Table II gives the results of the FEP/MD calculations with a complete breakdown of the various contributions to the total free energy. It is observed that the dominant contribution arises from the repulsive component from the aromatic side chain (−3.4 kcal/mol). As shown in Fig. 5 and in Table II, the repulsive component is actually larger when the selectivity filter is in its open-conductive state compared with the inactivated state. As a result, the ΔΔG for the F103A mutant favors the conductive state because of the larger repulsive free energy interaction of the side chain with the surroundings.

Thus, the final FEP/MD results are dominated by the unfavorable repulsive free energy of the active state relative to the inactive state:

where the thermodynamic coupling parameter s serves to switch on (s = 1) or off (s = 0) the repulsive contribution. The largest contribution to the difference in repulsive contributions arises when the difference in the derivative, $[∂ΔGrepF103A/∂s−∂ΔGrepWT(s)/∂s],$ reaches a maximum. Analysis of the configurations from FEP/MD trajectories shows that this occurs when the repulsive thermodynamic coupling parameter, s, is near 0.5. To clarify the origin of these contributions, it is useful to rewrite the thermodynamic derivatives in terms of explicit WCA-repulsive LJ interactions,

where $uαγrep(r;s)$ is the s-dependent soft-core WCA repulsive potential, and $ραγ(oo)(r;s)$ and $ραγ(oi)(r;s)$ are the average radial distribution functions between atom α (on Phe103) and another atom γ with the thermodynamic coupling parameter s for the open-conductive (oo) and inactivated states (oi), respectively. The WCA-repulsive LJ potential is a very short-range function. For this reason, the derivative $∂uαγrep(r;s)/∂s$ is nonzero only at distances that are shorter than or comparable to van der Waals contact. It follows that the integrand of Eq. 4 can contribute to ΔΔG_rep only if $ραγ(oo)(r;s)$ is larger than $ραγ(oi)(r;s)$ zero at very short distances. In other words, more atoms in the structure must overlap sterically with the aromatic side chain of Phe103 in the open-conductive (oo) state than in the inactivated (oi) state. This information is accessible only through this type of analysis; it cannot be deduced from the properties of the endpoint states (s = 0 and s = 1). Searching the configurations in the FEP/MD trajectories generated with s = 0.5 reveals that the dominant contributions arise from the interaction of Phe103 with Thr75 within the same subunit and Ile100 in the neighboring subunit (Fig. 1, B and C). Less important steric contacts include Leu105 in the same subunit and Val97 in the neighboring subunit. Thus, this analysis of FEP/MD trajectories supports previous conclusions based on the x-ray structure and equilibrium MD simulations (Cuello et al., 2010a).

A direct prediction of the analysis of FEP/MD trajectories is that a possible way to decrease the impact of Phe103 on the rate of inactivation might be to substitute a smaller residue at the position of Ile100. Interestingly, macroscopic pH-dependent currents for the KcsA-I100A mutant indicate that inactivation is indeed partially relieved in this channel. The experimentally observed effect of mutating an isoleucine at position 100 to an alanine on macroscopic pH-dependent currents is shown in Fig. 6. The macroscopic current response upon a pH jump shows that the I100A mutant inactivates more slowly and to a lesser extent than KcsA-WT. The mean inactivation time is ∼2.1 s for the WT channel and 8.5 s for the I100A mutant (not depicted). Although there is some variability (Fig. 6, red and blue traces), the slowdown in inactivation is accompanied by an increase in steady-state open probability to ∼30–40% for the mutant compared with 3–5% for KcsA-WT. The I100A mutant inactivates more slowly than WT on average, although there is some variability in value of the time constant. These observations are consistent with the Phe103–Ile100 intersubunit interaction noted in the analysis of the FEP/MD calculations. Although Ile100 is important, this experimental result confirms that Phe103 is the critical residue controlling the allosteric coupling between the filter and the intracellular gate. Although only the most unfavorable Phe103–Ile100 interaction can be relieved in the I100A mutation, it is likely that additional (smaller) unfavorable interactions made by Phe103 with other residues remain present. The latter are all removed in the F103A mutation, which may partly explain why the impact of the I100A mutation appears to be smaller than that of the F103A mutation.

Activation and inactivation gating are allosterically coupled. The closed conformation of the intracellular gate stabilizes the conductive conformation of the filter, whereas an open gate favors the inactivated filter. From the viewpoint of microscopic reversibility, however, any coupling of two structural elements must act in both directions, implying the possibility of “cross talk” between the intracellular gate and the selectivity filter. Thus, stabilizing the selectivity filter in its conductive state also promotes the closing of the intracellular gate. In KcsA, cross talk between the activation gate and the selectivity filter was first documented by electron paramagnetic resonance spectroscopy (Perozo et al., 1999), where opening of the inner gate leads to subtle conformational changes in the pore helix and outer vestibule. These conformational rearrangements, now known to relate to C-type inactivation at the selectivity filter, are also coupled to the opening of the activation gate (Chakrapani et al., 2007a,b). The present study has elucidated some of the molecular determinants responsible for coupling the state of the intracellular activation gate to the nonconductive inactivated conformation of the selectivity filter. Although only a single site-directed mutation was analyzed in detail using FEP/MD simulations, the results are encouraging and suggest that these types of computations could be a powerful approach for shedding light on the complex networks of residue–residue interactions underlying slow inactivation in K⁺ channels that have been revealed by previous experimental studies (Panyi and Deutsch, 2006).

Here, we have been mainly dealing with understanding the microscopic interactions involved with the thermodynamic stability/instability of the filter in the closed and open states. The mechanism by which the filter actually transitions from conductive to collapsed, which presumably would require a specification of the kinetic rate and the pathway, is beyond the present study. Identifying clearly the factors that affect the endpoints of such a path, however, is a necessary first step in elucidating the mechanism of inactivation. In particular, the present FEP/MD simulations are in relatively good agreement with experiments, which increases our confidence in the use of this type of calculations to study channel inactivation. Future work is underway to clarify the mechanism of recovery from inactivation using computational methods for determining the pathways connecting the different functional states of the channel (Pan et al., 2008; Gan et al., 2009).

This work was supported by the National Institute of Health through grant R01-GM62342 (to B. Roux), grant R01-GM57846 (to E. Perozo), and an NRSA postdoctoral fellowship (to A.C. Pan), and by the American Heart Association through grant 11SDG5440003 (to L.G. Cuello). The computations relied on the resources from the Extreme Science and Engineering Discovery Environment (XSEDE), which is supported by National Science Foundation (grant no. OCI-1053575).

Kenton J. Swartz served as editor.