Ion channels are membrane-bound enzymes whose catalytic sites are ion-conducting pores that open and close (gate) in response to specific environmental stimuli. Ion channels are important contributors to cell signaling and homeostasis. Our current understanding of gating is the product of 60 plus years of voltage-clamp recording augmented by intervention in the form of environmental, chemical, and mutational perturbations. The need for good phenomenological models of gating has evolved in parallel with the sophistication of experimental technique. The goal of modeling is to develop realistic schemes that not only describe data, but also accurately reflect mechanisms of action. This review covers three areas that have contributed to the understanding of ion channels: traditional Eyring kinetic theory, molecular dynamics analysis, and statistical thermodynamics. Although the primary emphasis is on voltage-dependent channels, the methods discussed here are easily generalized to other stimuli and could be applied to any ion channel and indeed any macromolecule.

Ion channels are membrane-bound enzymes whose catalytic sites are ion-conducting pores that open and close (gate) in response to specific environmental stimuli (voltage, ligand concentration, membrane tension, temperature, etc.; Hille, 2001). These stimuli activate highly specialized regulatory domains (an example is the voltage-sensing domain [VSD] found in many cation-selective channels) that are coupled to the pore gate through as yet incompletely understood mechanisms. Ion channels are important contributors to cell signaling and homeostasis and are strictly necessary for electrical conduction in nerve and muscle tissue. Our current understanding of channel gating is the product of over 60 years of voltage-clamp recording augmented by experimental intervention in the form of environmental, chemical, and mutational perturbations. Macroscopic ionic or capacitive gating currents reflecting the collective behavior of many channels can be interpreted with the help of descriptive (phenomenological) models that range from the standard two-state Boltzmann fit to more sophisticated analysis using multistate kinetic schemes (Hodgkin and Huxley, 1952; Vandenberg and Bezanilla, 1991; Zagotta et al., 1994; Schoppa and Sigworth, 1998b; Horrigan and Aldrich, 2002). Elementary gating events are resolvable using fluctuation analysis (Neher and Stevens, 1977; Sigworth, 1980, 1981; Conti and Stühmer, 1989; Sigg et al., 1994) or by measuring single-channel ionic currents (Colquhoun and Hawkes, 1995a). More recently, fluorometry and other optical techniques have been added to the repertoire of electrophysiological recordings (Mannuzzu et al., 1996; Cha and Bezanilla, 1997; Perozo et al., 1999). The tradition of modeling electrophysiological data dates back to the early 50s, when Hodgkin and Huxley published their celebrated 1952 papers culminating with a general description of the squid giant axon action potential (Hodgkin and Huxley, 1952). Alan Hodgkin writes about this time period (Hodgkin, 1976): “As soon as we began to think about molecular mechanisms it became clear that the molecular data would by itself yield only very general information about the class of system likely to be involved. So we settled for the more pedestrian aim of finding a simple set of mathematical equations which might plausibly represent the movement of electrically charged gating particles.”

What began as a quest to understand the physics of axon excitability (Hodgkin and Huxley had originally favored a carrier scheme) evolved into the slightly less ambitious objective of computing the time course of the action potential using mathematically defined units of activation, each containing four independent gating particles, thus anticipating the four-subunit structural motif now known to describe most ion channels. A plan to improve the fitting of voltage-clamp records by adding more gating particles was apparently tempered by the laborious necessity of performing numerical integration with a hand-cranked calculator: “Better agreement might have been obtained with a fifth or sixth power, but the improvement was not considered to be worth the additional complication” (Hodgkin and Huxley, 1952).

The past 60 years have seen tremendous advances in technique, including the ability to clone, mutate, and insert channels of any species into heterologous expression systems and the determination of x-ray crystal structures (Doyle et al., 1998; Long et al., 2005) in selected channels (and faster computers!). This has put us on the path to a residue-level understanding of ion channel structure and function. The need for good phenomenological models has evolved commensurately. More than two centuries after Luigi Galvani (1737–1798) formulated his “electrical hypothesis” of muscle contraction, we are approaching the threshold of a molecular and mechanistic description of membrane excitability. The goal is to develop gating schemes that do more than describe data: they also accurately reflect structures and mechanisms of action. This is an opportune time, then, to reflect on various kinetic and equilibrium models that have been used in the past, note recent advances, and comment on possible directions for the future of ion channel modeling.

The aim of this review is to explore common ground in the phenomenology of three fields that have contributed to the study of ion channels: traditional Eyring kinetic theory, statistical thermodynamics, and, most recently, molecular dynamics (MD) analysis. Much of the text will be devoted to the application of these fields as expressed by the following quantities: the transition rate constant (α), the potential of mean force (W [PMF]), and the partition function (Z). The emphasis will be on activated (barrier) processes with slow (millisecond) relaxation rates (readers wishing instead to review low-barrier electrodiffusion may start with these references: Bowman and Baglioni, 1984; Im and Roux, 2002). An essential component of the stated goal is a shared terminology. In view of the central role played by the voltage clamp in ion channel experimentation, the discussion of kinetics will use electrical analogues to the usual mechanical variables; for example, charge (q) ↔ distance, voltage (V) ↔ force, inductance (L) ↔ mass, and so on. However, recognizing that ion channels demonstrate tremendous variety in their responses to external stimuli, a subsequent section on thermodynamics is developed around canonical displacement-force pairs (ξ, Ω), thereby broadening the discussion from voltage-dependent (q, V) channels to ligand-gated (n, μ) and other types of channels. The scope of ion channel physiology is vast, and this review cannot hope to cover every facet of its phenomenology. As a former experimentalist turned itinerant modeler, the author does not claim a practicing knowledge in many areas touched upon here; in these an experimental viewpoint is adopted. Hopefully, the reader will be able to fill in any deficiencies from references spanning 80 years of theoretical and experimental work.

### Ion channel kinetics

Kinetic modeling of ion channels began with the Hodgkin–Huxley (HH) scheme (Hodgkin and Huxley, 1952) and continues unabated today with the widespread use of multistate Markov models and single-channel analysis. The “pedestrian aim” expressed by Hodgkin and Huxley, namely to simulate the shape and threshold nature of stimulated action potentials through simple kinetic equations, succeeded brilliantly, setting the stage for countless papers applying the HH formulism to model excitability in a variety of tissues. The building block of the HH model is the gating particle, which transitions between resting and activated states in a voltage-dependent manner. Because ion channels are composed of modular domains, each with its own specialization and distinct evolutionary origin (Schreiber et al., 1999; Anderson and Greenberg, 2001; Murata et al., 2005; Sasaki et al., 2006), the gating particle concept remains a viable one, although gating schemes have grown in complexity. The opening of a potassium (K+) channel in the original HH scheme required four randomly fluctuating “n” particles to be simultaneously active (the sodium channel description was similar but used three activating “m” particles and one inactivating “h” particle). K+ conduction in the HH model is therefore proportional to n4, with particle kinetics satisfying $n.=λ(n∞−n).$ The gating parameters λ = α + β and n = α/λ are functions of the voltage-dependent unidirectional rate constants α and β. The HH phenomenological equations are consistent with the notion that gating particles undergo transitions between macrostates (Fig. 1 A). The invention of the patch clamp (Neher and Sakmann, 1976) allowed single-channel experiments to confirm the existence of stochastic transitions between discrete conductance levels as predicted by the HH model.

The physics of particle transitions concerns the study of chemical escape from a metastable state, known as reaction rate theory, which has a long and interesting history (reviewed by Hänggi and Borkovec, 1990). The particular flavor of reaction rate theory best known to the ion channel community is the transition state theory (TST) of Eyring, which is premised on the somewhat tenuous assertion that reactants (R) rapidly thermalize with their surroundings until they reach the separatrix of the transition barrier (b), whereupon they inexorably turn to product (P). A one-dimensional bistable landscape serving as a model of unimolecular gating particle kinetics is shown in Fig. 1 A. The TST rate constant for the forward reaction is

$kTST=kThZb′ZRexp(−ΔE*kT),$

where the “universal prefactor” kT/h equals the reactive frequency of barrier crossing, the activation energy is ΔE* = EbER, and the Zs are statistical weights (partition functions) for the translational, vibrational, and rotational degrees of freedom of the reactant and barrier states (Eyring, 1935a,b). Pressure-volume work added to ΔE* yields the activation enthalpy ΔH* = ΔE* + Δv*P. The barrier partition function Zb is primed to indicate that it is missing the unstable forward translational component Ztran = (2πmkT)1/2/h per unit length that, when multiplied by the mean forward velocity (kT/2πm)1/2, is the source of the kT/h prefactor (Eyring, 1935b). The stable degree of freedom in ZR matching Ztran is the vibrational function Zvib = 2πkT/hωR, where ωR is the oscillator frequency of the reactant state. Separating Zvib from ZR reestablishes term balance in the ratio of statistical weights, yielding the more meaningful prefactor (ωR/2π)(Zb/ZR). The rate constant is further simplified by moving the Z ratio into the exponent, resulting in

$kTST=ωR2πexp(−ΔG*kT),$

where ΔG* = ΔH* − ΔS*T is the Gibbs activation energy, whose entropy ΔS* = kln(Zb/ZR) contains the nonreactive degrees of freedom. Contributions from external forces such as voltage V or chemical potential μ add new terms to ΔG* and generate the expanded free energy change ΔW* = ΔG* − Δq*V − Δn*μ −…, where Δq* (charge) and Δn* (binding number) are the canonical activation displacements to V and μ and other force-displacement pairs may be added.

The major failing of the TST (apart from Eyring’s misleading emphasis on the universal prefactor kT/h, a mere accounting trick) is that it does not properly account for the influence of the solvent or heat bath. Solvent interactions contribute not only additional degrees of freedom to ΔS*, but also generate fluctuations in the reacting coordinate that might reduce the rate of reaction (because not all barrier crossings lead to product). Eyring recognized the latter as a deficiency in his theory and suggested an empirical fix by multiplying kTST by a “transmission coefficient” κ. With these considerations, Eyring’s rate constant for the one-dimensional barrier assumes its final form:

$α=κkTST=κωR2πexp(−ΔW*kT).$
(1)

Despite its theoretical shortcomings, Eyring’s theory has been the mainstay of rate constant phenomenology in the ion channel literature. Voltage-dependent rate constants are usually expressed as variants of the following Eyring-like equations:

$α=αoexp(ΔqαVkT),$
(2a)
$β=βoexp(−ΔqβVkT),$
(2b)

where α and β are the forward (usually defined as positive charge carrying) and backward rates, respectively. The net charge displacement for the transition is Δq = Δqα + Δqβ. The preexponential factors in Eq. 2 (a and b) are related to Eyring’s theory through αo = κ(ωR/2π)exp(−ΔG*/kT) and βo = κ(ωP/2π)exp[−(ΔG* − ΔG)/kT], where ΔG = GPGR.

Channel gating is typically modeled by a network of discrete states, with each state connected to its neighbors through unimolecular rate constants. Such a network may be referred to as a discrete-state Markov (DSM) model. The word “Markovian” describes the absence of “memory,” in which predictions of future events depend solely on the present state. The number of states in DSM gating schemes is typically small (the HH model of the K+ channel consists of a linear arrangement of n = 5 states, and many empirical schemes do not exceed 10 states, although with allosteric models n may rise considerably). The existence of discrete closed-state transitions is supported by fluctuation analysis of gating currents (Conti and Stühmer, 1989; Sigg et al., 1994). Alternative gating models (Millhauser et al., 1988; Liebovitch, 1989) characterized by a “rough” energy landscape designed to simulate power-law kinetics seen in low-temperature protein motions (Frauenfelder et al., 1991) have not gained traction, primarily because the multiexponential kinetics predicted by DSM models adequately describe most electrophysiological data (McManus et al., 1988). The kinetics of an n-state DSM model are obtained by solving the master equation dp/dt = pA, where p(t) is the 1 × n vector of state probabilities and A is an n × n array of rate constants αij known as the “Q-matrix” (Colquhoun and Hawkes, 1995b). Knowing A, one can derive the nonstationary mean value and flux for any observable a through the vector relations 〈a〉 = 〈p|a〉 and $〈a.〉$ = 〈pA|a〉, where a is the n × 1 state vector of a values, and p(t) is solved as an initial value problem through numerical integration or eigenvector decomposition. The latter method yields

$p=∑r=0n−1〈p(0)|vr〉urexp(−λrt),$

where ur and vr are the left and right eigenvectors of A and the corresponding eigenvalues −λr are the macroscopic decay rates. In a closed (gating) system satisfying detailed balance, all λr are distinct, real, and positive, except for the zeroth eigenvalue λ0 = 0, whose left eigenvector u0 is the equilibrium probability distribution p(∞). Eigenvalues for the simple HH K+ channel scheme are integer multiples of λ = α + β. In contrast, an open system (pore, transporter) may decay with damped oscillations (Re[λr] ≥ 0) toward a nonequilibrium steady-state (Schnakenberg, 1976). An important element of Q-matrix theory is kinetic analysis of single-channel data, but this subject is too vast to cover here (reviewed by Qin, 2007).

Although DSM models have been immensely helpful in studying gating at the network level, the physical origin of the rate constant, specifically the basis for the Eyring variables κ, Δq, ω, and ΔG*, remains poorly understood. One obstacle to experimentally determining rate constants is the difficulty in measuring single transition rates without the aid of empirical models. An important step forward was achieved when crystal structures of ion channels were made available (Doyle et al., 1998; Long et al., 2005), providing the starting point for all-atom MD simulations and the theoretical ability to compute rate constants in silico. Although most of the MD work on ion channels has focused on ion permeation in pores (Roux and Karplus, 1991; Allen et al., 1999; Shrivastava and Sansom, 2000; Bernèche and Roux, 2001), long-time scale trajectories of voltage sensor activation have recently been published (Delemotte et al., 2011; Jensen et al., 2012), setting the stage for detailed MD studies of gating.

Because of the large time scale difference between the step size (approximately femtoseconds) of MD integration and intervals between ion hops in pores (nanoseconds to microseconds) or the relaxation time of a gating particle (approximately milliseconds), special techniques are required to compute rate constants from MD trajectories. The progenitor of these techniques is the reactive flux method of Chandler (1978), which utilizes as the source of calculation an ensemble of short MD trajectories originating at the barrier peak. With the right choice of reaction coordinate q and provided the precise location of the q barrier is known, the reactive flux of the forward reaction kf(t) can be computed from

$kf(t)=〈q.(0)δ[q(0)−qb]H[q(t)]〉〈1−H(q)〉,$
(3)

where brackets signify stationary averaging. The delta function δ[q(0) − qb] restricts the numerator to trajectories starting at the barrier peak at time t = 0, and the Heaviside function H[q(t)] (=1 if q(t) > qb, otherwise =0) ensures that only trajectories on the product side are counted. The denominator 〈1 − H(q)〉 is the equilibrium probability of the reactant state. At time t = 0+, one obtains kf = kTST because only positive velocities can enter the product side. However, after a brief (approximately picoseconds) time period τmol during which one or more barrier crossings may occur, the fate of the trajectory is determined (R or P), and the reactive flux averaged over many trajectories plateaus to a value equal to the phenomenological rate constant α (Eq. 1). Thus, the ratio of reactive flux taken at t = 0+ and t = τmol yields the transmission coefficient κ = α/kTST. It is worth nothing that neither κ nor kTST is independent of the choice of reaction coordinate or the location of the transition state, only the product α = κkTST holds that distinction, but a poor choice of reaction coordinate leads to poor statistics and possibly an erroneous estimate of α (Berne et al., 1988; Crouzy et al., 1994). Subsequent developments in reactive flux theory have improved computational speed (White et al., 2000) and broadened the application to energy landscapes with complicated topographies lacking a well-defined transition state (Dellago et al., 1999).

Although successfully predicting the value of a rate constant by MD simulation is a notable achievement, a more informative picture of the transition process is obtained by mapping out the energy conformational landscape with phenomenological models. Similar to how experimentalists interpret voltage-clamp data using DSM models, molecular dynamicists can project the complete phase space (system plus heat bath) sampled by MD simulation onto a small number of degrees of freedom described by a low-dimensional stochastic differential equation. The generalized Langevin equation (GLE) is such an equation. Given a reaction coordinate q with mass L, the GLE describes the motion of q influenced by three forces: the gradient of a free energy landscape W(q), also known as the PMF, a dissipating force involving a time-dependent friction coefficient R(q,t), and a rapidly fluctuating force Λ(t). The one-dimensional GLE has the following form (Straub et al., 1987):

$Lq..=−W′(q)−R(q,t)*q.+Λ(q,t),$

where dots and primes denote time and spatial derivatives, respectively. The convolution operation

$R*q.=∫0tR(q,s)q.(t−s)ds$

“filters” the velocity $q.$ through the friction kernel R, introducing memory into the short time-scale dynamics. Λ is generally assumed to be Gaussian and is therefore completely specified by the first and second moments: 〈Λ〉 = 0 and 〈Λ(0)Λ(t)〉 = kTR(t). The latter equation establishes a relation between the friction-damping and fluctuating forces known as the second fluctuation-dissipation theorem (Kubo, 1966).

The task of designing MD protocols to compute W(q) and R(q,t) specific to electrophysiological processes has been systematically undertaken over the past two decades, much of it through the efforts of the Roux laboratory (Crouzy et al., 1994; Roux, 1997, 1999; Bernèche and Roux, 2001, 2003; Im and Roux, 2001). Recent work has sought to improve upon traditional MD force fields by incorporating polarization effects that play an important role in determining the strength of electrostatic interactions (Lopes et al., 2009; Bucher and Rothlisberger, 2010). Regardless of MD methodology, the process of coarse-graining trajectory data onto the GLE begins with the identification of a reaction coordinate q. The obvious (though not necessarily the most efficient for sampling purposes) choice of q in voltage-dependent transitions is the gating (or ionic) charge displacement, which hereafter will be understood to be synonymous with the reaction coordinate q. Determining the charge displacement as a function of configuration requires computing the detailed membrane potential profile in the channel with fixed (i.e., nonmobile and nonpolarizable) charges turned off (Roux, 1997). The PMF, correct to first order in V, is then given by W(q) = G(q) − qV (Roux, 1999, 2008). The zero-voltage PMF G(q) is obtained from G(q) = −kTlnp(q), where p(q) is the probability histogram of q acquired over a long MD trajectory at constant simulated temperature and pressure. Sampling statistics for the infrequently visited barrier region are improved by confining the simulation to sequential regions (umbrella windows) along the activation landscape (Torrie and Valleau, 1977; Roux, 1995). The same umbrella-sampled simulations can be used to compute R(q,t) from the stationary velocity autocorrelation function $〈q.(0)q.(t)〉〈q.2〉−1$ (Straub et al., 1987).

Extracting the GLE parameters W(q) and R(q,t) from MD simulations is not a trivial undertaking, particularly if the reaction pathway (which may be multidimensional) is not well known. This is certainly true for the case of gating transitions in ion channels because available crystal structures provide only frozen “snap shots” of a particular configuration, and mechanisms of gating are poorly understood. For this reason, the early focus has been on ion permeation, where the inconvenience of dealing with an open system (grand canonical ensemble) is compensated by the existence of a relatively straightforward reaction pathway, where permeant ions have well-defined masses (Roux and Karplus, 1991). To compute rate constants from the landscape parameters W(q) and R(q,t), dynamic Langevin simulations can be performed (Lange and Grubmüller, 2006; Gordon et al., 2009), or the Eyring parameters ωR, ΔW*, and κ can be extracted, allowing α to be calculated from Eq. 1. To execute the latter, one generally resorts to approximating the energy landscape by a piecewise sequence of parabolic basins and barriers (Fig. 1 A), treating R as constant. Each segment (basin or barrier) of the landscape behaves as an oscillator satisfying the simplified GLE equation:

$q..+γ*q.+εω2q=ξ(t),$
(4)

where the memory kernel γ(t) = R(t)/L and the frequencies ω of local oscillations are equal to (LC)−1/2. The fluctuating force ξ(t) = Λ(t)/L satisfies 〈ξ(0)ξ(t)〉 = (kT/L)γ(t). The inverse “spring constant” C is related to the curvature of the basin or barrier through C−1 = εW″(q), where ε is 1 for (stable) basins and −1 for (unstable) barriers. Damped oscillations in the reactant and product basins (ε = 1) are exactly described by an electrical RLC circuit, where R is resistance, L is inductance, and C is capacitance. It is straightforward to demonstrate that for a linear voltage drop over a distance D, the variables R, L, C, and q are related to their mechanical equivalents (friction, mass, spring constant, and distance) through factors of z/D, where z is the charge valence of the oscillating degree of freedom.

With Eq. 4 governing the reactant (εR = 1) and barrier (εb = −1) segments of the PMF, we easily derive ωR = (LCR)−1/2 and ΔW* = ΔG* − Δq*V. Obtaining the transmission coefficient κ from Eq. 4 (known as the linearized GLE because the restoring force is linear in q) is more involved, but several derivations have been published (Grote and Hynes, 1980; Hänggi and Mojtabai, 1982; Pollak, 1986), including an elegant treatment based on reactive flux described by Eq. 3 (Kohen and Tannor, 2000). The satisfyingly simple result, attributed to Grote–Hynes (GH), is

$κ=λωb.$
(5)

In Eq. 5, λ is the largest positive root of $ϑ−1=s2+γ^s−ωb2,$ where ϑ(s) is the transfer function of Eq. 4 within the barrier region and

$γ^(s)=∫0∞γ(t)exp(−st)dt$

is the Laplace transform of γ(t). Although evaluating λ requires specific knowledge of γ(t), it is instructive to examine limiting cases. When coupling forces to the bath are weak, implying a small value of γ(t) or the combination of a narrow barrier and/or small L leads to a large ωb, then the rapidly accelerating q at the barrier is decoupled from slower frictional forces, and the positive root of ϑ−1 is λ = ωb, leading to the TST result κ = 1. (In reality the TST value represents an upper limit to the rate constant because a transition from spatial and velocity diffusion to energy diffusion occurs as R is reduced to very low values, which again diminishes the value of κ [Kramers, 1940; Mel’nikov and Meshkov, 1986], although the so-called small friction regime [Fig. 1 B] is irrelevant under physiological conditions.) In contrast, the combination of a broad barrier, large L, and strong friction coupling (all probable characteristics of a gating process) reduces the relative contribution of the ωb2 term in ϑ−1, leading to a smaller value of λ and enabling $γ^(s=λ)$ to approach its zero-frequency value

$γ^(0)=∫0∞γ(t)dt.$

Under these so-called intermediate-to-large friction conditions, the positive root of ϑ−1 is easily obtained:

$λ=γ^(0)2(1+(2ωbγ^(0))2−1).$
(6)

If $γ^(0)>>ωb,$Eq. 6 reduces to $λ=ωb2/γ^(0),$ leading to the large friction (LF; Smoluchowski) limit $κ=ωb/γ^(0)$ and the following well-known expression for the rate constant (Kramers, 1940):

$α=ωRωb2πγ^(0)exp(−ΔW*kT).$
(7)

Kramers actually derived Eqs. 6 and 7 by implicitly assuming rapid (“memory less”) friction damping, which we can invoke by specifying γ(t) = 2γδ(t). In this case $γ^=γ$ for all frequencies, reducing Eq. 4 to the conventional Langevin equation:

$q..+γq.+εω2q=ξ(t),$
(8)

where 〈ξ(0)ξ(t)〉 = 2kT(γ/L)δ(t). Substituting the now time-independent γ for $γ^(0)$ in Eqs. 7 and 8, we obtain the explicit forms of Kramers’ intermediate-large (Kr) and LF expressions for α (plotted as a fraction of kTST in Fig. 1 B). Contraction to a pure spatial diffusion (classical Brownian motion) occurs in the LF limit where the system is always at terminal velocity, enabling one to neglect the acceleration term in Eq. 8. This yields, after dividing by γ,

$q.+εqτ=F(t),$
(9)

where τ = γ/ω2 = RC is a novel relaxation time that was implicated in the early fast component of Shaker K+ channel gating currents (Sigg et al., 2003). The fluctuating force in Eq. 9 satisfies 〈F(0)F(t)〉 = 2Dδ(t), where the diffusion constant D is related to the friction R through D = kT/R (known as the Einstein relation). Evaluating Eq. 9 for the barrier region at constant R again yields $λ=ωb2/γ,$ from which the following forward and backward rate constants can be derived, in analogy to Eq. 2 (a and b) and consistent with Eq. 7:

$α=12πRCRCbexp(−ΔW*kT),$
(10a)
$β=12πRCPCbexp(−(ΔW*−ΔW)kT).$
(10b)

An expression analogous to Eq. 10b was used to evaluate the reverse transition in the first barrier of the Shaker K+ channel activation sequence (Sigg et al., 2003).

The oscillator dynamics described by Eq. 9 can be generalized to arbitrary landscape topologies W(q), which is described by the LF Smoluchowski equation, shown here in Fokker–Planck form (Van Kampen, 1992):

$∂p∂t=∂∂q(W′(q)R(q)p+kTR(q)∂p∂q),$
(11)

where p(q,t) is a probability distribution and R(q) is allowed to vary across q. The behavior of Eq. 11 applied to various landscape models of gating has been extensively described in Sigg et al. (1999). Operating in the spatial Smoluchowski regimen has advantages over the use of “ballistic” equations (such as the GLE), which feature velocity as well as spatial diffusion and thus require a reduced mass L. The absence of mass in Eqs. 9 and 11 benefits gating models because it may be difficult to define the masses of regulatory domains functioning as gating particles. Provided the energy barrier in an arbitrary {W(q), R(q)} Smoluchowski landscape is sufficiently large (>5 kT), then the forward transition rate constant is equal to the inverse of the mean first passage time (Van Kampen, 1992), which is easily computed from the following double integral (refer to Fig. 1 A for limits of integration):

$α−1=∫qRqPdqR(q)kTexp(W(q)kT)∫q(−)qdq′exp(−W(q′)kT).$
(12)

Direct evidence for the validity of the LF assumption in gating dynamics is still lacking, although in one study of primitive gating in the gramicidin channel there was good agreement between Kramers’ LF theory and experimental rates (Crouzy et al., 1994). There is more experience with MD studies of barrier hopping of permeant and blocking ions in pores. The results from three such studies are summarized in Table 1. In each case, Langevin dynamics were found to be overdamped with κ ranging from 0.1 to 0.25. In two studies (Roux and Karplus, 1991; Tolokh et al., 2002), comparisons of κ values were made between the “gold-standard” reactive flux method (Eq. 3), the GH equation (Eq. 5), and Kramers’ intermediate-to-large friction formula (Eq. 7). There was excellent agreement between the reactive flux and GH methods, whereas Kramers’ theory underestimated κ by 20–40%, suggesting memory friction plays a nonnegligible role in ion permeation. In contrast, it is likely, based on our earlier arguments, that the larger scale of conformational changes associated with gating implies Markovian dynamics, and in the presence of extensive external (solvent, lipid) and internal (protein) damping, the Kramers–Smoluchowski condition (Eqs. 7 and 912) should apply. In an interesting side note, comparing Na+ conduction in the gramicidin pore (Roux and Karplus, 1991) with Ba2+ block in the KcsA channel (Rowley and Roux, 2013), the computed 107-fold difference in kTST was primarily caused by differences in barrier height ΔG* because in both cases the reactant state frequency ωR/2π was within an order of magnitude (!) of the universal rate constant kT/h ≈ 6 × 1012 s−1.

The preceding discussion sought to bridge macroscopic (DSM and TST) and microscopic (GLE, Langevin, and Smoluchowski) equations of motion by focusing on an area of common ground: the rate constant value as determined by the Eyring, reactive flux, GH, and Kramers formulas. Confidence in the validity of a phenomenological description is increased if an experimentally determined rate is consistent with the prediction of MD. However, kinetic modeling may not necessarily be the most efficient method of assigning a gating mechanism to experimental findings. To take a specific example, the mechanism by which VSDs gate the pore in voltage-dependent ion channels remains unclear. Barring direct visualization of the gating process (a distinct future possibility given advances in MD!), the researcher must infer the nature of VSD–pore coupling using indirect methods, for example by fitting ionic or gating currents to a DSM model (Vandenberg and Bezanilla, 1991; Schoppa and Sigworth, 1998b; Rothberg and Magleby, 2000) or by identifying key residues of interaction (typically highly conserved across channels) through mutation analysis or other means of functional disruption (Yang and Horn, 1995; Holmgren et al., 1998; Horn et al., 2000; Tiwari-Woodruff et al., 2000; Arcisio-Miranda et al., 2010). Model-independent methods such as mutant cycle analysis (Horovitz, 1996; Yifrach and MacKinnon, 2002; Zandany et al., 2008) can quantitatively measure coupling between residues of interacting domains such as the pore and the VSD provided that the free energies associated with the observed process are accurately measured. Mutant cycle analysis relies on conservation of energy and displacement and the ability to sum these extensive variables across a series of perturbations. Relaxation rates (eigenvalues) obtained from a multistate ion channel are not additive, although small perturbations in the activation energy of a single rate constant do add, a fact which has been used to an advantage in studying mutant cycles (Schreiber and Fersht, 1995). The problem is that, as previously stated, determining rate constants from single-channel analysis is resource intensive and typically model dependent (Colquhoun and Hawkes, 1995a; Qin, 2007). Extreme conditions are necessary to measure elementary rates directly from macroscopic currents (Sigg et al., 2003; Chakrapani and Auerbach, 2005). As a result, it is common to analyze perturbation effects through the use of standard equilibrium curves of observable quantities such as conductance (G) and gating charge (Q). This falls under the scope of thermodynamics, and although DSM models have equilibrium conditions baked into their structure, there exist simpler and more powerful methods of analyzing equilibrium curves that yield the desired energies and displacements, often in model-independent fashion. This is the subject of the next two sections.

### Thermodynamic models

In modeling electrophysiological data, thermodynamics has traditionally taken a backseat to kinetic theory, which through the use of DSM schemes enjoys widespread familiarity, but is clumsy at handling thermodynamic derivations. The thermodynamic counterpart to the Q-matrix is the partition function Z, a statistically weighted sum (or integral) of states from which all macroscopic quantities are derived (discussed in any statistical mechanics textbook; e.g., Hill, 1960). For a one-dimensional PMF, the partition function is

$Z=∫exp(−W(q)/kT)dq,$

where the integral is taken over all q. A normalizing factor for Z is unnecessary because it merely adds a constant value to the principle quantity of interest, the channel’s chemical potential Φ, obtained from the “bridge” equation Φ = −kTlnZ. The Gibbs–Duhem equation, which describes the interdependence of intensive quantities in a thermodynamic system (Waldram, 1985), relates dΦ to incremental changes in system constraints (T, P, V, μ,…, Ω) through: dΦ = −〈sdT + 〈vdP − 〈qdV − 〈ndμ −…− 〈ξ〉dΩ, where canonical displacements 〈s〉, 〈v〉, 〈q〉, 〈n〉,…, 〈ξ〉 are written in lower case to indicate division by N (channels). In this way, we can derive the value of each displacement using 〈ξ〉 = −∂Φ/∂Ω = kT∂lnZ/∂Ω (an exception of sorts is entropy, which is not typically thought of as a “displacement” and possesses an extra term: 〈s〉 = kT∂lnZ/∂T + klnZ). For example, the mean gating charge displacement is obtained from

$〈q〉≡kT∂lnZ/∂V=Z−1∫qexp(−W(q)/kT)dq,$

the derivation of which requires energy to be linearly related to charge: W(q) = G(q) − qV. We’ll examine the thermodynamics of a piecewise-harmonic two-state landscape in some detail before discussing more complex systems.

One of the many convenient properties of Z is the unrestricted ability to subpartition it into smaller parts, which we’ll refer to as the parsing principle. In a conformational energy landscape with large barriers, we can divide the landscape into stable (positive curvature) and unstable (negative curvature) regions. For example, the partition function of the bistable potential in Fig. 1 A may be written Z = ZR + Zb + ZP, where

$ZX=∫{X}exp(−W(q)/kT)dq$

is the subpartition function for the region {X} bounded by turning points in the curvature G″(q). The local potentials are second-order expansions around local extrema (qX, GX) with G″ = ±CX−1 and additionally subject to a linear voltage term −qV. A large barrier minimizes the contribution of Zb, effectively yielding ZZR + ZP. After integrating and choosing (qR, GR) as the zero reference, we obtain

$ZR=2πkTCRexp(CRV22kT),$
(13a)
$ZP=2πkTCPexp(CPV22kT)exp(−ΔWkT).$
(13b)

In Eq. 13b, ΔW = ΔG − ΔqV, where, as previously, we define ΔG = GPGR and Δq = qPqR.

We also find it useful to propose a modified partition function that emphasizes the barrier region:

$Zb*=∫{b}exp(W(q)/kT)dq,$

which differs from the usual definition of Z by assigning a positive sign to the exponent. Performing the integral, we obtain:

$Zb*=2πkTCbexp(CbV22kT)exp(ΔW*kT),$
(14)

where ΔW* = ΔG* − Δq*V.

The preceding exercises lay the groundwork for a rational method of discretizing a continuous energy landscape into a finite set of basin and barrier coordinates (qX, ΦX). These coordinates are derived from the relations: qX = ±kT∂lnZX/∂V and ΦX = ∓kTlnZX, where the upper sign applies to basins and the lower sign to barriers. Applying these expressions to Eqs. 13 and 14 yielded the coordinates of the discretized landscapes in Fig. 1 A (dashed lines). Several results ensue from this coarse-graining procedure. Revisiting kinetic theory, we see that in evaluating Eq. 12 for the bistable potential at constant R, the inner integral is almost complete and approximates ZR at the value of q where the outer integral begins to increase in value. After factoring out ZR, the outer integral is roughly equal to Zb*. We therefore conclude that

$α≈(kTR)1ZRZb*=Dexp(−ΔΦ*kT),$
(15)

where ΔΦ* = Φb − ΦR and the Einstein relation was again used. Inserting Eqs. 13a and 14 into Eq. 15 and neglecting V2 terms that introduce small capacity corrections to the barrier height, we see that the result is consistent with Eq. 10a obtained from Kramers’ LF theory. Also, Δq* formally equals Δqα from Eq. 2a. The right side of Eq. 15 has a pleasing simplicity. The preexponential factor is the diffusion constant D = kT/R, and the activation energy ΔΦ* is easily determined for any shape of the barrier landscape. Detailed examination of the fast component of gating current in Shaker yielded an estimate of D = 41 eo2/ms, corresponding to the early stages of activation (Sigg et al., 2003).

Returning to the subject of thermodynamics, we derive the average charge displacement for the bistable potential using 〈q〉 = pRqR + pPqP, where pX = ZX/Z are state probabilities of the stable R and P regions. Evaluating the response to a change in voltage ΔV with the aid of Eq. 13 (a and b), we obtain two components of charge movement (Fig. 1 A):

$〈q〉=ZR(V)CR+ZP(V)CPZR(V)+ZP(V)ΔV+ZP(V)ZR(V)+ZP(V)Δq.$
(16)

The first term in Eq. 16 is a small capacity charge δq = pRδqR + pPδqP, where the “drift” charges δqX = CXΔV would be experimentally indistinguishable from the linear membrane capacitance if not for the fact that basin states vary in width (Sigg et al., 2003). We recall from earlier that the “fast” gating current described by δq decays with time constant RC.

The second term in Eq. 16 relates to the much slower (relaxation time λ−1) “hopping” or transition charge Δq, which is easily distinguished experimentally from the fast component, both from a temporal standpoint (λ−1 >> RC) and also in magnitude: Δq >> |δqP − δqR| (the minus sign in the fast charge |δqP − δqR| is a peculiarity of the linear subtraction methods used in gating current measurements). Thus, in describing 〈q〉, we typically neglect the fast charge in favor of the larger transition charge, whose equilibrium distribution is 〈q〉 = [K/(1 + K)]Δq, where KpP(∞)/pR(∞) = ZP/ZR is the equilibrium constant. Substituting Eq. 13 (a and b) into the above expression for K, again neglecting the V2 terms corresponding to the fast capacity charge, we obtain

$K=CPCRexp(−ΔWkT),$
(17)

which agrees with Eq. 10 (a and b) and the detailed balance requirement: αpR(∞) = βpP(∞).

We note the parallel between 〈q〉/Δq = K/(1 + K) and the venerable Boltzmann function B = KB/(1 + KB), where KB = exp(ψB/kT) and the bias potential ψB = ξB(Ω − ΩB) is a function of the force-displacement Boltzmann parameters ΩB and ξB. The practice of Boltzmann fitting is both widely used and often abused, as shoehorning the complexity of a multidomain channel into two states requires rigorous justification (for example, Yifrach and MacKinnon, 2002). Nevertheless, Boltzmann fitting is a convenient bookkeeping technique for tracking changes to activation curves in response to systemic perturbations. Nonsymmetrical activation curves can be fitted with m Boltzmann functions, although strictly speaking the physical basis underlying such a procedure is the existence of m independently gating domains described by Z = (1 + K1)f1(1 + K2)f2…(1 + Km)fm, where the fi values are population fractions. Faced with choosing between Boltzmann fitting or analysis with an n-state Markov model, the Boltzmann method may in many cases be the more sensible alternative (see commentary by Shem-Ad and Yifrach, 2013). Nevertheless, a systematic approach that (a) measures thermodynamically relevant quantities, (b) uses fewer adjustable variables than a comparable kinetic model, and (c) acknowledges known or proposed structural and functional relationships between channel components is a welcomed addition to the modeling armamentarium. Such an approach (linkage analysis) has recently been applied to channel gating (Chowdhury and Chanda, 2010, 2012a, 2013; Sigg, 2013), but before discussing it, we should review ion channel thermodynamics a bit more.

One consequence of the parsing principle is the various ways Z may be expressed for a multistate channel. Four parsing methods (mathematically equivalent to alternative factoring schemes) applied to a linear three-state model are shown in Fig. 2. These fall into three general categories based on the method of computing the equilibrium mean for an arbitrary observable a:

$absolute: 〈a〉=∑i=1npiai,pi=∂lnZ∂lnZi=ZiZ;$
(18a)
$relative: 〈a〉=a1+∑i=2nθiΔai,θi=∂lnZ∂lnΥi,Yi=ZiZ1,Δai=ai−a1;$
(18b)
$relational: 〈a〉=a1+∑KϕKΔaK,ϕK=∂lnZ∂lnK,K=exp(−ηKkT).$
(18c)

The absolute and relative methods (Eq. 18, a and b) sum over n and n − 1 states, respectively, and are useful in deriving thermodynamic relations (Conti, 1986; Sigg and Bezanilla, 1997). The relational method (Eq. 18c) differs from the other two in that the sum runs over individual transitions rather than over all states. For a linear scheme, this is no great advantage because there exist up to n − 1 distinct transitions in an n-state model, but for allosteric models possessing a small number of transitions, the dimensionality of the system can be substantially reduced. For example, the classic Monod–Wyman–Changeux (MWC) model (Monod et al., 1965) has 10 states but only two unique transitions (and one allosteric factor). A 70-state model of the BK channel (Horrigan and Aldrich, 2002) is described by three transitions and three allosteric factors. The quantity ϕK in relational parsing is the equilibrium curve for the transition K with range {0...nK}, where nK is the number of identical K particles. Channel opening is typically described by a pore gate with equilibrium constant L coupled to other processes (J, K). The open probability Po in such cases is equal to ϕL = ∂lnZ/∂lnL, which varies from 0 to 1. It is generally straightforward to compute ϕK for any transition K because for DSM models, Z can always be expanded as a polynomial of K.

The “particle potential” ηK in Eq. 18c, defined as the potential difference ΔΦK between reactant and product states of a particle, deserves special attention because it contains local thermodynamic information related to K activation. Contributions to ηK may include particle changes in volume (ΔvK) and entropy (ΔsK), charge movement (ΔqK), stoichiometries (ΔnK) of chemical reactions (e.g., ion, ligand, or solvent binding), and, relevant to allosterically regulated channels, energies of interaction WKJ between K and other particles J. Written explicitly, ηK = ΔEK − ΔsKT + ΔvKP − ΔqKV − ΣΔnKμK +…+ ΣWKJ. If one wishes to emphasize a particular extensive parameter ΔξK, one can factor it out of the sum and regroup the remaining terms in the form of a corresponding half-activation value ΩK, yielding ηK = −ΔξK(Ω − ΩK). For example, rewriting Eq. 17 for the equilibrium constant of a two-state system in a way that emphasizes gating charge, one obtains: K = exp(−ηK/kT) = exp[Δq(VVK)/kT], where VK = Δq−1G − (kT/2)ln(CP/CR)]. If positive Ω activates a particle (the usual convention), then we can substitute ψK for −ηK, where ψK = ΔξK(Ω − ΩK) is the bias potential introduced earlier in the context of the Boltzmann function.

The last of the four categories in Fig. 2 is the hierarchical approach, which parses Z in order of most to least relevant equilibrium constants (L is the obvious top choice given the central role of the pore). Hierarchical parsing can be seen as a relational method that leaves open the possibility of deeper categorization as new functionalities are elicited (Sigg, 2013).

In light of the exponential dependence of an equilibrium constant K to its bias potential ψK = ΔξK(Ω − ΩK), a plot of lnA versus Ω, in which A = Na〉 is a marker of K activation, should yield information about Ω-linked processes relevant to K, specifically the canonical displacement ΔξK and coupling energies WKJ linking K to other Ω-dependent particles J. This is the basis for site-specific (local) linkage analysis developed by Jeffries Wyman to study cooperative oxygen binding in hemoglobin (Wyman, 1967) and recently adapted to ion channels (Chowdhury and Chanda, 2010, 2012a, 2013; Sigg, 2013). The parameter of interest in local linkage analysis is the Hill energy, operationally defined for an observable A as WH[a] = kTln[(AAmin)/(AmaxA)]. Defining the probability of K activation by PKPA = (AAmin)/(AmaxAmin) and rearranging, we see that WH[a] = kTln[PA/(1 − PA)] is the change in negative system energy Φ from K activation: −ΔΦK = kTlnPAkTln(1 − PA). Provided the particle of interest K is weakly coupled to other gating particles, a plot of WH[a] versus Ω is sigmoidal with linear lower and upper asymptotes that run parallel to ψK but are vertically separated by the total coupling energy WKΩ between K and the various J particles. An example of a local variable is the pore conductance G. In principle, a plot of WH[g] versus V should yield model-independent estimates of the intrinsic pore charge ΔqL and the total energy of pore–VSD coupling (Sigg, 2013).

Wyman also described “global” linkage, which concerns energies of interaction between mediators of dissimilar forces (Ω = V, T, u1, u2, etc.; Wyman, 1967; Chowdhury and Chanda, 2013; Sigg, 2013). The parameter of interest here is the total energy of activating Ω-sensitive particles: ΔΦΩ = ΔξmaxΩM, where Δξmax is total displacement related to Ω, and ΩM, which has been called the median (or mean, depending on interpretation) force of activation, and is defined as the value of Ω that divides the plot of 〈ξ〉 versus Ω into equal areas (the median interpretation; Wyman, 1967; Chowdhury and Chanda, 2012a), or equivalently, as the average Ω with respect to capacitance d〈ξ〉/dΩ (the mean interpretation; Di Cera and Chen, 1993; Sigg, 2013). The voltage parameter VM was extensively discussed in Chowdhury and Chanda (2012a), where they demonstrated that for multistate channels marked differences may exist between VM and the Boltzmann-fit parameter VB.

The model-independent displacements (ΔξK, Δξmax) and energies (ΔΦK, ΔΦΩ) obtained from local and global linkage plots are characteristic thermodynamic quantities and, when possible, should be obtained in lieu of corresponding Boltzmann parameters ξB and ξBΩB (Chowdhury and Chanda, 2012a; Bezanilla and Villalba-Galea, 2013; Sigg, 2013). The Boltzmann equation has been used as an all-purpose tool to fit any activity curve (most typically G-V or Q-V), but its validity not does extend beyond n > 2. Linkage analysis is theoretically preferable where an accurate measurement of energy or displacement is critical to the interpretation of results. A pertinent example introduced earlier is the method of double mutant cycle analysis, in which the nonadditivity of single mutant perturbation energies in double or higher dimensional mutants is considered strong evidence of interaction between mutated residues (Carter et al., 1984; Horovitz, 1987). Amino acid residues separated by 6 Å or less have been shown to interact with energies ranging from 0.5 to 7 kcal/mol (Schreiber and Fersht, 1995), with electrostatic interactions contributing to larger energies. Applied to ion channels, double mutant cycle analysis with Boltzmann fitting has been used to measure pore–VSD interactions (Yifrach and MacKinnon, 2002; Zandany et al., 2008). Chowdhury and Chanda (2010) have suggested that a local (ΔΦK based) linkage analysis, which they call χ analysis, be considered as a universally applicable method for quantifying mutational effects, although to date no practical application of the technique has been published, at least in part because of the technical difficulty of measuring logarithmic-scale sensitivities of local markers of activation. Chowdhury and Chanda have more recently presented a global (ΔΦΩ based) method of quantifying pore–VSD interactions based on the Q-V curve that was used to quantify energetic interactions in a triad of interacting residues in the Shaker K+ channel (Chowdhury, S., et al. 2014. Biophysical Society 58th Annual Meeting. 3748-Pos Board B476).

Before turning to allosteric channels, a subject ideally suited to linkage analysis, we highlight a key result from Conti (1986) that has languished in relative obscurity within the ion channel literature for nearly 30 years but covers much of the thermodynamic theory described here. Starting from the equilibrium relations

$〈a〉=∑ipiai,$

pi = Zi/Z, and ξi = kT∂lnZi/∂Ω, Conti derived

$kT∂ln〈a〉∂Ω=cov(a,ξ)〈a〉,$
(19)

where cov(a,ξ) = 〈aξ〉 − 〈a〉〈ξ〉 is the covariance between any state-dependent observable a and the Ω-linked displacement ξ. Eq. 19 was derived for discrete states but also covers diffusion landscapes as a consequence of the parsing principle. In the special case of a = ξ, Eq. 19 simplifies to the well-known relation for capacitance: kT∂〈ξ〉/∂Ω = 〈Δξ2〉. The Conti relation is insensitive to the number of channels N because replacing 〈a〉 with A = Na〉 leaves Eq. 19 unchanged. We can equate the open probability Po with the unit-normalized G-V curve = (GGmin)/(GmaxGmin), where G = Ng〉 and Gmin is a baseline conductance. Applying Eq. 19 to Po with Ω = V, we obtain an expression for the “activation charge” 〈qa〉 ≡ kTdln(GGmin)/dV = 〈qo〉 − 〈q〉, where 〈qo〉 = 〈Poq〉/Po is the mean charge displacement among open states (Sigg and Bezanilla, 1997). Similarly, defining the closed probability Pc = 1 − Po, we obtain the complementary quantity 〈qd〉 ≡ −kTdln(GmaxG)/dV = 〈q〉 − 〈qc〉, where 〈qc〉 = 〈Pcq〉/Pc is the closed-state charge displacement, and we also have 〈q〉 = Pcqc〉 + Poqo〉. Finally, we recognize that for positive voltage-activating channels the slope of the conductance Hill plot m = dWH[g]/dV ≡ 〈qa〉 + 〈qd〉 = 〈qo〉 − 〈qc〉 asymptotically approaches ΔqL for extreme negative (m(−) ≈ 〈qa〉) and positive (m(+) ≈ 〈qd〉) voltages because under these limiting conditions, the voltage sensor is locked into either the resting or activated position and the pore can gate quasi-independently. Channels for which m(−) = m(+) measured under experimental conditions may be said to exhibit “weak” pore–VSD coupling, and gating schemes can either conform to this condition or violate it through “strong” or “obligatory” coupling, as discussed presently.

### Allosteric models of gating

Some 10 years after Hodgkin and Huxley’s groundbreaking work on the action potential, a different, nonkinetic approach yielding the well-known MWC allosteric model was applied to the problem of oxygen binding in hemoglobin. The MWC scheme proposed that the O2 affinities of the n = 4 hemoglobin subunits increased simultaneously with the flip of a two-state L parameter (Monod et al., 1965). The ligand-gated ion channel community was quick to catch on, with publication two years later of an n = 2 MWC model for the acetylcholine channel (Karlin, 1967). The physical basis of L in both hemoglobin and ligand-gated channels is thought to lie with symmetry-preserving residue interactions at the interfaces between subunits (Changeux and Edelstein, 2005). In voltage-dependent ion channels, a case can be made equating the pore with the “L particle” because pore opening is a binary process n-fold connected to neighboring VSDs. It wasn’t until the 1990s that allosteric models of voltage-gated channels were considered seriously, culminating in two influential papers from the Aldrich laboratory. The first (Ledwell and Aldrich, 1999) studied a late opening transition in the voltage-sensitive Shaker K+ channel by mutating three noncharged residues of the mid-to-distal portion of the VSD S4 segment that interact with the pore domain (Pathak et al., 2007). The authors found that this so-called ILT mutant (Smith-Maxwell et al., 1998) led to single-exponential time courses in the ionic current and a far right-shifted G-V curve, consistent with a rate-limiting opening transition. They interpreted their findings in the context of traditional obligatory gating models of Shaker activation (Bezanilla et al., 1994; Zagotta et al., 1994; Schoppa and Sigworth, 1998b) in which pore opening occurs only at the end of the activation sequence, as originally prescribed by the HH scheme. The apparent absence of more than one open state in obligatory gating is at odds with the notion of weak allosterism, but a kind of strong coupling is inferred by the requirement that all four VSDs must activate before pore opening. Support for obligatory gating in Shaker and other homologous voltage-gated K+ channels came from “limiting slope” experiments (Almers, 1978; Sigg and Bezanilla, 1997), in which the activation charge 〈qa〉 experimentally approaches the value of total activation charge Δqmax at very low Po ∼10−7 (Islas and Sigworth, 1999). At such low Po, the gating charge versus voltage curve (Q-V) is negatively saturated (〈q〉 ≈ 0) and 〈qa〉 = m(−) = Δqmax). Evidently, the open state charge 〈qo〉 remains maximal in Shaker even at very low open probabilities, as demonstrated by 〈qo〉 ≡ 〈qa〉 + 〈q〉 ≈ Δqmax. This rules out the existence of intermediate-charged open states, for which 〈qo〉 < Δqmax. Obligatory gating in Shaker is consistent with single-exponential open time distributions in singles analysis (Hoshi et al., 1994; Schoppa and Sigworth, 1998a) and the observation that just one charged voltage sensor is sufficient to gate the pore shut (Gagnon and Bezanilla, 2009). However, if one takes the position that obligatory gating represents the strong coupling limit of an underlying allosteric scheme (a position which is admittedly at odds with current conventional wisdom), then one obtains a lower estimate of the pore–VSD coupling energy W required to approach the limiting slope m(−) = Δqmax for Po = 10−7 by assuming that the corresponding point in WH[g] occurs halfway between lower and upper asymptotes (corresponding to the steepest portion of the Hill plot), yielding W = 2kTln(10−7) ≈ 32 kT, or ∼8 kT per subunit. An allosteric model of Shaker consistent with the ILT experiments that nearly achieves this coupling requirement is presented later.

The second of the aforementioned Aldrich papers (Horrigan and Aldrich, 2002), which focused on the large conductance (BK) Ca2+- and voltage-activated K+ channel, is arguably the most complete characterization of allosterism in a voltage-gated cation channel to date. The energies of interaction between the VSD (J particle), the pore (L particle), and the Ca2+-binding domain (K particle) obtained in this landmark study ranged from 0.5 to 1.9 kcal/mol (0.9 to 3.2 kT) per subunit, with the upper limit in energy corresponding to J-L coupling. The negative limiting slope m(−) in BK for Po < 10−3 (high [Ca2+]) and 10−7(low [Ca2+]) yielded ΔqL = 0.3 eo. Weak allosterism has been implicated in gating models of other channels, with derived interaction energies ranging from 1.0 to 1.2 kcal/mol (1.7 to 2.0 kT) in the voltage- and temperature-gated TRPM8 channel (Brauchi et al., 2004), 0.9 to 1.8 kcal/mol (1.6 to 3.2 kT) in voltage- and cyclic nucleotide-gated HCN channels (Altomare et al., 2001; Chen et al., 2007; Ryu and Yellen, 2012), 0.4 to 1.0 kcal/mol (0.7 to 1.8 kT) in the voltage-gated KCNQ1 channel (Osteen et al., 2012), 2.9 kcal/mol (5.0 kT) in a sodium channel (Arcisio-Miranda et al., 2010), and 1.6 kcal/mol (2.7 kT) in a calcium channel (Marks and Jones, 1992).

The basic unit of interaction in an allosteric scheme is the pairwise coupling between two gating particles that we’ll refer to as a “cooperon.” Consider an isolated cooperon between J and L particles whose Gibbs energy landscapes are shown in Fig. 3 A. The corresponding single particle partition functions are ZJ = 1 + Jo and ZL = 1 + Lo, where Jo and Lo represent native equilibrium constants. Interactions between the activation states of J and L, if not too strong or distorting, generate a two-dimensional landscape of four basins arranged in a cycle (Sigg and Bezanilla, 2003). The partition function in each basin state can be constructed from the particle states (j, l) and their strengths of interaction Djl. Thus, the cooperon partition function is Z = D00 + JoD10 + LoD01 + JoLoD11. The interaction free energies and charges are given by Wjl = −kTlnDjl and δqjl = kT∂lnDjl/∂V, with the latter usually set to zero. Dividing Z by D00 and reassigning terms yields the more compact expression Z = 1 + J + L + JLD, where J = JoD10/D00, L = LoD01/D00, and D = D00D11/D10D01 are renormalized variables derived from native equilibrium constants and raw coupling factors (Chowdhury and Chanda, 2010). There are no restrictions on manipulations of this type because, as noted earlier, normalizing factors do not change energy relations between states. Because the global coupling factor D targets the double product state A-O (Fig. 3), we refer to the renormalized scheme as “product coupled.” The fact that D is a composite of all four native interactions implies that an apparent absence of interaction between two particles (D = 1) does not necessarily imply complete particle independence, but may be a consequence of balanced energies, specifically WW00 + W11W10W01 = 0. We also consider the possibility of excess charge δq = kT∂lnD/∂V, which renders D voltage dependent provided that δq = δq00 + δq11 − δq10 − δq01 ≠ 0. Evidence has been offered in support of a small (0.2 eo) δq in internal coupling between BK voltage sensor particles (Pantazis et al., 2010).

In most voltage-gated ion channels, Po increases with voltage. A physical mechanism consistent with this property is when voltage sensors assist the opening of a reluctant pore gate through positive interaction (negative W) in the A-O state. This “pull-open” mechanism is illustrated on the right side of Fig. 3 B; it satisfies the product-coupled partition function Z = 1 + J + L + JLD, with D > 1. However, one can also imagine a “push-closed” mechanism (Fig. 3 B, left) in which the pore opens intrinsically, but only after activation of the J particle relieves a steric hindrance D′ in the R-O state, which prevents pore opening while J is at rest (Chowdhury and Chanda, 2012b; Horrigan, 2012). The push-closed mechanism is described by Z = 1 + J + LD′ + JL′ with D′ < 1. It is entirely equivalent to the pull-open mechanism provided we set D′ = 1/D and L′ = LD. The conductance Hill plots of the two mechanisms coincide if their respective L particle biases ψL and ψL′ are vertically separated by the coupling energy W = −kTlnD (Fig. 3 B). In either case, if we allow L to become very small and D very large, then Z ≈ 1 + J + JL′ (Arcisio-Miranda et al., 2010), which describes the three-state obligatory model of Fig. 2 and represents the limiting case of L-J coupling in Fig. 3 B (dotted red line). The extension to a multi-cooperon scheme consisting of a central pore particle L radially coupled to n J particles is straightforward, where Z = (1 + J)n + L(1 + JD)n for the pull-open mechanism, Z = (1 + J)n + L′(D′ + J)n with L′ = LDn for the push-closed mechanism, and Z = (1 + J)n + LJn for the (n + 2)-state obligatory model.

If, practically speaking, one could eliminate the source of cooperativity D (or D′) from a product-coupled scheme, then the ensuing position of VL (or VL′) on the V axis should favor a particular gating mechanism (pull-open or push-closed). It is in this context that an alternative interpretation of the ILT experiments is proposed, one that regards Shaker gating as the manifestation of strong coupling rather than following a strictly obligatory scheme. We entertain the possibility that the ILT mutation completely uncouples the pore from the voltage sensor, reducing the value of D to the neutral value 1.0. In fact, some residual pore–VSD coupling may exist in the ILT mutant (Smith-Maxwell et al., 1998), but we’ll treat this as negligible. The observed G-V shift to the right in ILT experiments should therefore favor the pull-open hypothesis. The difficulty with this interpretation relates to the change in the Q-V curve in response to ILT, which undergoes an opposite shift to the left (Ledwell and Aldrich, 1999). It is generally the case that changes to a single coupling factor or equilibrium constant is inconsistent with opposing activation shifts (Chowdhury and Chanda, 2012b). Similar examples of diverging G-V and Q-V curves have since been observed in a handful of Shaker mutants, mostly involving residues found in the interaction region between the VSD and the pore (segments S5 and S6 in the pore and the S4–S5 linker; Soler-Llavina et al., 2006; Haddad and Blunck, 2011). Interestingly, neutralizing salt bridge–forming charged residues in the S2 and S3 segments of the VSD (Papazian et al., 1995; Tiwari-Woodruff et al., 1997) have also generated opposing shifts (Seoh et al., 1996).

Muroi et al. (2010) provided a plausible mechanism to explain the divergence of G-V and Q-V curves upon observing the phenomena in sodium channel mutants targeting the domain III VSD–pore linkage region. They showed that if a mutation were to simultaneously affect separate coupling factors in the L-J cooperon, specifically if one factor targets the doubly “reactant” state (R-C) and the other the “product” state (A-O), a scheme referred to here as “doubly coupled,” then relieving both sources of coupling could generate opposing shifts in ϕL and ϕJ and also therefore in the G-V and Q-V curves. The concept of interactions favoring “agreement” states such as R-C and A-O is a common one in statistical mechanics. It explains the quaternary symmetry requirements of the MWC model (Monod et al., 1965) and underlies the Lenz–Ising model of lattice interactions (Niss, 2005).

Incorporating double coupling into an n = 4 radial cooperon scheme (Fig. 4) yields the partition function

$Z=(C+J)4+L(1+JD)4,$
(20)

where, as before, L and J are the pore and voltage sensor equilibrium constants, and the coupling factors associated with R-C and A-O configurations are C = exp(−WC/kT) and D = exp(−WD/kT), respectively. Assuming in Shaker a total charge per channel Δqmax = 13 eo (Aggarwal and MacKinnon, 1996; Noceti et al., 1996; Seoh et al., 1996; Islas and Sigworth, 1999), Eq. 20 adequately describes the changes made by the ILT mutations to the G-V and Q-V curves using just five adjustable parameters (ΔqL, VL, VJ, WC, and WD), all related to displacements and energies obtained from linkage analysis (Fig. 4, B and C).

Although the n = 4 doubly coupled scheme is allosteric, with half of its 32 states conducting, the ILT fits account for obligatory-like gating in Shaker (Fig. 4 C) by distributing the total energy of interaction Wtot = 4|WC + WD| ≈ 26 kT equally among the four subunits, with 70% of Wtot allotted to R-C interactions (4|WC| ≈ 18 kT) and 30% to A-O interactions (4|WD| ≈ 8 kT). This value of Wtot is roughly 80% of the earlier back-of-the-envelope estimate of 32 kT for the minimum VSD–pore coupling energy required to emulate obligatory gating in Shaker, leaving open the possibility that a more completely uncoupled mutant (see Haddad and Blunck, 2011) might yield a larger value of Wtot through greater separation of G-V and Q-V curves. We note that mathematically there is no difference between stabilizing the double reactant and product configurations and destabilizing intermediate states. Therefore, ILT or similarly acting mutations could conceivably act by weakening favorable interactions in agreement states R-C and A-O as suggested above or by relieving “strain” in intermediate states R-A and C-O (Pathak et al., 2005; Haddad and Blunck, 2011).

One could argue that, by factoring out the quantity C4 from Eq. 20 to obtain Z = (1 + J′)4 + L′(1 + JCD)4, with effective coupling constant CD and redefined equilibrium constants J′ = J/C and L′ = L/C4, one simply reverts back to a product-coupled scheme of the pull-open or reformulated push-closed type. This is true, but the above function differs from those of the earlier schemes in that both J′ and L′ are decreased in value with increasing C (recall that J in the pull-open or push-closed schemes was unchanged in the limit of large coupling). Hence, if C → ∞, the doubly coupled scheme reduces to a two-state model described by Z = 1 + LJ4(D/C)4, provided we also have a large D to maintain energy balance (note the absence of opening bias if D = C). A limiting doubly coupled scheme of this type might explain the binary activation of the four-subunit pore gate, which prefers all S6 helical bundle segments to be in the closed or open state, thereby satisfying quaternary symmetry required by MWC theory (Chapman and VanDongen, 2005).

A mechanistic interpretation compatible with the doubly coupled scheme is if disruption of packing interactions between the pore and VSD is necessary for intermediate gating transitions to occur, a conjecture supported by the findings of long MD trajectories of K+ channel activation (Jensen et al., 2012). Although the doubly coupled scheme is overly simplistic, neglecting for example the multistep activation of the VSD (Zagotta et al., 1994; Schoppa and Sigworth, 1998b; Gamal El-Din et al., 2010; Delemotte et al., 2011; Haddad and Blunck, 2011; Yarov-Yarovoy et al., 2012), it adopts a thermodynamic viewpoint in which the parameters of interest consist solely of energies and displacements. This lines up perfectly with the way we often think of mechanisms of action, which focus on mechanical details such as which residues in domain X interact with those in domain Y and how these interactions are modified when there is a change in conformation. Kinetic modeling provides additional information about the time course of events, but this may not be of central importance in determining structure–function relationships, and the added complexity of a kinetic model can obscure interpretation. The double coupling framework is easily expanded to include multistep VSD transitions and could be used to explore other features of gating, such as intersubunit pore domain interactions and channel inactivation. Adding new gating particles and transitions to a multimeric allosteric model can result in an unwieldy proliferation of kinetic states, although equilibrium curves can still easily be calculated using the partition function formulism. Of course, a successful model must be able to predict gating and ionic currents, and for this purpose it is necessary to adopt a systematic approach to kinetic analysis of allosteric models, as discussed next.

### Kinetics in allosteric models

Constructing a kinetic model of allosterism can be a daunting task. Thermodynamics offers only incomplete guidance on constraints placed by allosteric interactions on the perturbation of rate constants, and the number of kinetic states rises geometrically with cooperon number. A reasonable solution to the first problem is the concept of linear free energy relationships (LFERs; Leffler, 1953; Hammond, 1955; Grosman et al., 2000; Azaria et al., 2010; Edelstein and Changeux, 2010). The basic idea is that if, structurally speaking, the transition state (b) is intermediate between the reactant (R) and product (P) states, then the perturbation of state potentials ΦX by an arbitrary force Ω should satisfy the weighted average: ∂Φb/∂Ω = x∂ΦP/∂Ω + (1 − x)∂ΦR/∂Ω, where the weighting factor x (LFER parameter) is a number between 0 and 1. Rearranging terms, we have: ∂ΔΦ*/∂Ω = x∂ΔΦ/∂Ω, where ΔΦ* = Φb − ΦR and ΔΦ = ΦP − ΦR. We obtained earlier (Eq. 15) α = Dexp(−ΔΦ*/kT), from which we derive, assuming Φ is linear in Ω, α(Ω)/α(0) ≡ exp[−(∂ΔΦ*/∂Ω)Ω/kT] = exp[−x(∂ΔΦ/∂Ω)Ω/kT]. For the familiar case of Ω = V and substituting ∂ΔΦ/∂V = −Δq, we obtain α = αoexp(xΔqV/kT), where, comparing with Eq. 2a, we establish that Δqa (=Δq*) = xΔq. Extending the above argument to both forward and backward rates, we can write $α=αoΓVx$ and $β=βoΓVx−1$, where ΓV = exp(ΔqV/kT) is the voltage perturbation factor. The ratio α/β yields the usual expression for the equilibrium constant K = Koexp(ΔqV/kT), where Ko = αoo. We can introduce additional factors for other applied forces (ΓT, ΓP, Γμ,…) as well as for allosteric interactions ($DJj,DLl,$…), where j, l,… are the number of activated J, L,… particles. Generally speaking, the LFER parameter x of a particular perturbation ∂Ω depends on the position of the transition peak along the coordinate ξ = −kT∂Φ/∂Ω. However, we can make the zero-order assumption that x has the same value for any perturbation and further specify a transition-specific relaxation rate ν = αo(1−x)βox, allowing us to economically write: α = νKx and β = νKx−1, where $K=KoΓTΓPΓVΓμ…DJjDLl....$ Expressing this as a particle potential ηK = −kTlnK, we have

$ηk=ΔEK−ΔsKT+ΔvKP−ΔqKV−ΔnKμ−…−jWKJ−lWKL−....$
(21)

Eq. 21 is a complete thermodynamic description of K activation that includes allosteric interactions with other particles. After incorporating Eq. 21 into the above expressions for K, α, and β, the resulting relations represent the basis for a phenomenological rate description. Combinatorial factors (σK = 1…n) must be included as part of a global kinetic scheme if there are n identical K particles. A heuristic model of the BK channel comprising m = 4 particle species arranged in fourfold symmetry around the pore for a total of 17 particles and 7,752 configurational states (twice the number of ways n = 4 subunits can occupy 2m = 16 subunit states) was recently evaluated in this manner (Sigg, 2013), although the large state space made kinetic analysis computationally intensive, requiring Monte Carlo methods.

A shortcut bypassing the rigorous Q-matrix procedure of kinetic analysis would be of considerable benefit when dealing with large-scale allosteric models. The earlier trick of computing equilibrium curves by summing over elementary transitions using the relational form of the partition function (Eq. 18c) suggests one might construct an analogous set of coupled HH-like kinetic equations, one for each gating particle, where the discrete coupling numbers j, l,… in Eq. 21 are replaced with the time-varying ϕJ, ϕL,… in the form of a mean-field approximation. Unfortunately, the Q-matrix of allosteric models is not so easily coarse grained. Although the mean-field procedure satisfies equilibrium requirements, it distorts the time course of relaxation even for the simple cooperon scheme. However, if one particle is rate limiting, a successful approach analogous to the Born–Oppenheimer strategy from quantum chemistry can be used. Consider an n-meric channel with a rate-limiting pore L relative to other regulatory particles (νL << νJ, νK,…). The partition function can be written as Z = Zcn + LZon, where the closed-state component Zc = f(J, K,…) is a polynomial function of subunit equilibrium constants J, K,…. The same polynomial describes the open state partition function Zo = f(JD, KC,…) except that the equilibrium constants are multiplied by allosteric factors D, C,…. The pore then activates in the manner of an HH particle, characterized by $ϕ.L=λL(ϕL∞−ϕL),$ where ϕL = ∂lnZ/∂lnL and λL = αL + βL. The rate constants αL and βL are Boltzmann-weighted sums of elementary pore rates across {J, K,…} substates, each individually constructed in the manner detailed in the previous paragraph. Evaluating the two sums yields the compact expression

$λL=ZL1−x(ZxZcZo)nνL,$
(22)

where Zx = f(JDx, KCx,…) is a sliding partition function which varies from Zc to Zo in the range x = {0…1}. The fast gating particles are “slaved” to L and therefore relax with rate λL after an initial rapid decay. An approach mathematically analogous to Eq. 22 was used to compute ionic current relaxation rates in the BK channel, whose pore kinetics are roughly 10 times slower than those of its regulatory domains (Cox et al., 1997; Horrigan and Aldrich, 2002).

### Conclusion

The history of modeling ion channels is filled with many successes and some dead ends. Even the dead ends are useful in that they test the limits of our understanding. This review has attempted to bridge microscopic and macroscopic views of ion channel dynamics through an exploration of common ground approaches to the phenomenology of gating. Although the lion’s share of attention was directed toward voltage-dependent channels, the methods discussed here are easily generalized to other stimuli and are applicable to any channel and indeed any macromolecule. Looking forward, close cooperation between the disciplines of MD and electrophysiology will likely play an important role in unraveling the molecular mechanism of gating. Along these lines, a consensus view of the activation pathway for voltage-gated ion channels based on decades of experimentation and simulation modeling has recently been published (Vargas et al., 2012). A logical next step would be to map out potentials of mean force for the gating landscape under physiological conditions. This would require a large effort (some might call it a pipe dream) but offers a substantial payoff: the validation, based on a channel’s molecular structure, of specific gating models used with confidence to interpret experimental data, with the ultimate aim of improved drug development and alleviation of disease. A combined approach has recently been applied (Ostmeyer et al., 2013) to the question of why recovery from C-type inactivation in the K+ pore is so slow (∼5–20 s). The answer? Three water molecules trapped behind the selectivity filter of each subunit bolster the pinched conformation of the inactivated pore. PMF calculations in this study demonstrated a large energy barrier (∼25 kcal/mol) to pore conduction unless all (or most) of the 12 waters were released into solution. With the benefit of hindsight and using thermodynamic arguments, one might have guessed the role of inactivating waters based on earlier work that studied the pressure dependence of slow inactivation in Shaker (Meyer and Heinemann, 1997), but molecular simulation combined with kinetic modeling was necessary to illustrate the precise mechanism in nearly cinematic detail. We should expect more such exciting developments arising from the interface between molecular and macroscopic dynamics. Whether an accurate phenomenology of gating dynamics follows the script suggested by the contents of this review, or whether new modeling paradigms come into play, based in part on advances in computer simulations but informed by experiment, should become clearer over the next decade or so as more light is shed on this fascinating and important issue.

Helpful comments on the manuscript by Dr. Baron Chanda and Dr. Richard Aldrich are gratefully acknowledged.

The author declares no competing financial interests.

Kenton J. Swartz served as editor.

Aggarwal
S.K.
,
MacKinnon
R.
.
1996
.
Contribution of the S4 segment to gating charge in the Shaker K+ channel
.
Neuron.
16
:
1169
1177
.
Allen
T.W.
,
Kuyucak
S.
,
Chung
S.H.
.
1999
.
Molecular dynamics study of the KcsA potassium channel
.
Biophys. J.
77
:
2502
2516
.
Almers
W.
1978
.
Gating currents and charge movements in excitable membranes
.
Rev. Physiol. Biochem. Pharmacol.
82
:
96
190
.
Altomare
C.
,
Bucchi
A.
,
Camatini
E.
,
Baruscotti
M.
,
Viscomi
C.
,
Moroni
A.
,
DiFrancesco
D.
.
2001
.
Integrated allosteric model of voltage gating of HCN channels
.
J. Gen. Physiol.
117
:
519
532
.
Anderson
P.A.
,
Greenberg
R.M.
.
2001
.
Phylogeny of ion channels: clues to structure and function
.
Comp. Biochem. Physiol. B Biochem. Mol. Biol.
129
:
17
28
.
Arcisio-Miranda
M.
,
Muroi
Y.
,
Chowdhury
S.
,
Chanda
B.
.
2010
.
Molecular mechanism of allosteric modification of voltage-dependent sodium channels by local anesthetics
.
J. Gen. Physiol.
136
:
541
554
.
Azaria
R.
,
Irit
O.
,
Ben-Abu
Y.
,
Yifrach
O.
.
2010
.
Probing the transition state of the allosteric pathway of the Shaker Kv channel pore by linear free-energy relations
.
J. Mol. Biol.
403
:
167
173
.
Berne
B.J.
,
Borkovec
M.
,
Straub
J.E.
.
1988
.
Classical and modern methods in reaction rate theory
.
J. Phys. Chem.
92
:
3711
3725
.
Bernèche
S.
,
Roux
B.
.
2001
.
Energetics of ion conduction through the K+ channel
.
Nature.
414
:
73
77
.
Bernèche
S.
,
Roux
B.
.
2003
.
A microscopic view of ion conduction through the K+ channel
.
100
:
8644
8648
.
Bezanilla
F.
,
Villalba-Galea
C.A.
.
2013
.
The gating charge should not be estimated by fitting a two-state model to a Q-V curve
.
J. Gen. Physiol.
142
:
575
578
.
Bezanilla
F.
,
Perozo
E.
,
Stefani
E.
.
1994
.
Gating of Shaker K+ channels: II. The components of gating currents and a model of channel activation
.
Biophys. J.
66
:
1011
1021
.
Bowman
C.L.
,
Baglioni
A.
.
1984
.
Application of the Goldman-Hodgkin-Katz current equation to membrane current-voltage data
.
J. Theor. Biol.
108
:
1
29
.
Brauchi
S.
,
Orio
P.
,
Latorre
R.
.
2004
.
Clues to understanding cold sensation: thermodynamics and electrophysiological analysis of the cold receptor TRPM8
.
101
:
15494
15499
.
Bucher
D.
,
Rothlisberger
U.
.
2010
.
Perspectives on: Molecular dynamics and computational methods: Molecular simulations of ion channels: a quantum chemist’s perspective
.
J. Gen. Physiol.
135
:
549
554
.
Carter
P.J.
,
Winter
G.
,
Wilkinson
A.J.
,
Fersht
A.R.
.
1984
.
The use of double mutants to detect structural changes in the active site of the tyrosyl-tRNA synthetase (Bacillus stearothermophilus)
.
Cell.
38
:
835
840
.
Cha
A.
,
Bezanilla
F.
.
1997
.
Characterizing voltage-dependent conformational changes in the Shaker K+ channel with fluorescence
.
Neuron.
19
:
1127
1140
.
Chakrapani
S.
,
Auerbach
A.
.
2005
.
A speed limit for conformational change of an allosteric membrane protein
.
102
:
87
92
.
Chandler
D.
1978
.
Statistical mechanics of isomerization dynamics in liquids and the transition state approximation
.
J. Chem. Phys.
68
:
2959
2970
.
Changeux
J.-P.
,
Edelstein
S.J.
.
2005
.
Allosteric mechanisms of signal transduction
.
Science.
308
:
1424
1428
.
Chapman
M.L.
,
VanDongen
A.M.J.
.
2005
.
K channel subconductance levels result from heteromeric pore conformations
.
J. Gen. Physiol.
126
:
87
103
.
Chen
S.
,
Wang
J.
,
Zhou
L.
,
George
M.S.
,
Siegelbaum
S.A.
.
2007
.
Voltage sensor movement and cAMP binding allosterically regulate an inherently voltage-independent closed-open transition in HCN channels
.
J. Gen. Physiol.
129
:
175
188
.
Chowdhury
S.
,
Chanda
B.
.
2010
.
Deconstructing thermodynamic parameters of a coupled system from site-specific observables
.
107
:
18856
18861
.
Chowdhury
S.
,
Chanda
B.
.
2012a
.
Estimating the voltage-dependent free energy change of ion channels using the median voltage for activation
.
J. Gen. Physiol.
139
:
3
17
.
Chowdhury
S.
,
Chanda
B.
.
2012b
.
Perspectives on: Conformational coupling in ion channels: Thermodynamics of electromechanical coupling in voltage-gated ion channels
.
J. Gen. Physiol.
140
:
613
623
.
Chowdhury
S.
,
Chanda
B.
.
2013
.
Free-energy relationships in ion channels activated by voltage and ligand
.
J. Gen. Physiol.
141
:
11
28
.
Colquhoun
D.
,
Hawkes
A.G.
,
.
1995a
.
The principles of the stochastic interpretation of ion-channel mechanisms
.
In
Single-Channel Recording.
Second edition.
Sakmann
B.
,
Neher
E.
,
editors
.
Springer
,
New York
.
397
482
.
Colquhoun
D.
,
Hawkes
A.G.
,
.
1995b
.
A Q-matrix cookbook: How to write only one program to calculate the single-channel and macroscopic predictions for any kinetic mechanism
.
In
Single-Channel Recording.
Second edition.
Sakmann
B.
,
Neher
E.
, editors.
Springer
,
New York
.
589
633
.
Conti
F.
,
1986
.
The relationship between electrophysiological data and thermodynamics of ion channel conformations
.
In
Ion channels in neural membranes: Proceedings of the 11th International Conference on Biological Membranes Held at Crans-sur-Sierre, Switzerland, June 10–14, 1985.
Ritchie
J.M.
,
Keynes
R.D.
,
Bolis
L.
,
editors
.
A.R. Liss
,
New York
.
25
41
.
Conti
F.
,
Stühmer
W.
.
1989
.
Quantal charge redistributions accompanying the structural transitions of sodium channels
.
Eur. Biophys. J.
17
:
53
59
.
Cox
D.H.
,
Cui
J.
,
Aldrich
R.W.
.
1997
.
Allosteric gating of a large conductance Ca-activated K+ channel
.
J. Gen. Physiol.
110
:
257
281
.
Crouzy
S.
,
Woolf
T.B.
,
Roux
B.
.
1994
.
A molecular dynamics study of gating in dioxolane-linked gramicidin A channels
.
Biophys. J.
67
:
1370
1386
.
Delemotte
L.
,
Tarek
M.
,
Klein
M.L.
,
Amaral
C.
,
Treptow
W.
.
2011
.
Intermediate states of the Kv1.2 voltage sensor from atomistic molecular dynamics simulations
.
108
:
6109
6114
.
Dellago
C.
,
Bolhuis
P.G.
,
Chandler
D.
.
1999
.
On the calculation of reaction rate constants in the transition path ensemble
.
J. Chem. Phys.
110
:
6617
6625
.
Di Cera
E.
,
Chen
Z.Q.
.
1993
.
The binding capacity is a probability density function
.
Biophys. J.
65
:
164
170
.
Doyle
D.A.
,
Morais Cabral
J.
,
Pfuetzner
R.A.
,
Kuo
A.
,
Gulbis
J.M.
,
Cohen
S.L.
,
Chait
B.T.
,
MacKinnon
R.
.
1998
.
The structure of the potassium channel: Molecular basis of K+ conduction and selectivity
.
Science.
280
:
69
77
.
Edelstein
S.J.
,
Changeux
J.-P.
.
2010
.
Relationships between structural dynamics and functional kinetics in oligomeric membrane receptors
.
Biophys. J.
98
:
2045
2052
.
Eyring
H.
1935a
.
The activated complex in chemical reactions
.
J. Chem. Phys.
3
:
107
115
.
Eyring
H.
1935b
.
The activated complex and the absolute rate of chemical reactions
.
Chem. Rev.
17
:
65
77
.
Frauenfelder
H.
,
Sligar
S.G.
,
Wolynes
P.G.
.
1991
.
The energy landscapes and motions of proteins
.
Science.
254
:
1598
1603
.
Gagnon
D.G.
,
Bezanilla
F.
.
2009
.
A single charged voltage sensor is capable of gating the Shaker K+ channel
.
J. Gen. Physiol.
133
:
467
483
.
Gamal El-Din
T.M.
,
Heldstab
H.
,
Lehmann
C.
,
Greeff
N.G.
.
2010
.
Double gaps along Shaker S4 demonstrate omega currents at three different closed states
.
Channels (Austin).
4
:
93
100
.
Gordon
D.
,
Krishnamurthy
V.
,
Chung
S.-H.
.
2009
.
Generalized Langevin models of molecular dynamics simulations with applications to ion channels
.
J. Chem. Phys.
131
:
134102
134112
.
Grosman
C.
,
Zhou
M.
,
Auerbach
A.
.
2000
.
Mapping the conformational wave of acetylcholine receptor channel gating
.
Nature.
403
:
773
776
.
Grote
R.F.
,
Hynes
J.T.
.
1980
.
The stable states picture of chemical reactions. II. Rate constants for condensed and gas phase reaction models
.
J. Chem. Phys.
73
:
2715
2732
.
G.A.
,
Blunck
R.
.
2011
.
Mode shift of the voltage sensors in Shaker K+ channels is caused by energetic coupling to the pore domain
.
J. Gen. Physiol.
137
:
455
472
.
Hammond
G.S.
1955
.
A correlation of reaction rates
.
J. Am. Chem. Soc.
77
:
334
338
.
Hänggi
P.
,
Borkovec
M.
.
1990
.
Reaction-rate theory: fifty years after Kramers
.
Rev. Mod. Phys.
62
:
251
341
.
Hänggi
P.
,
Mojtabai
F.
.
1982
.
Thermally activated escape rate in presence of long-time memory
.
Phys. Rev. A.
26
:
1168
1170
.
Hill
T.L.
1960
.
An Introduction to Statistical Thermodynamics
.
,
New York
.
508 pp
.
Hille
B.
2001
.
Ion Channels of Excitable Membranes
. Third edition.
Sinauer Associates
,
Sunderland, MA
.
814 pp
.
Hodgkin
A.L.
1976
.
Chance and design in electrophysiology: an informal account of certain experiments on nerve carried out between 1934 and 1952
.
J. Physiol.
263
:
1
21
.
Hodgkin
A.L.
,
Huxley
A.F.
.
1952
.
A quantitative description of membrane current and its application to conduction and excitation in nerve
.
J. Physiol.
117
:
500
544
.
Holmgren
M.
,
Shin
K.S.
,
Yellen
G.
.
1998
.
The activation gate of a voltage-gated K+ channel can be trapped in the open state by an intersubunit metal bridge
.
Neuron.
21
:
617
621
.
Horn
R.
,
Ding
S.
,
Gruber
H.J.
.
2000
.
Immobilizing the moving parts of voltage-gated ion channels
.
J. Gen. Physiol.
116
:
461
476
.
Horovitz
A.
1987
.
.
J. Mol. Biol.
196
:
733
735
.
Horovitz
A.
1996
.
Double-mutant cycles: a powerful tool for analyzing protein structure and function
.
Fold. Des.
1
:
R121
R126
.
Horrigan
F.T.
2012
.
Perspectives on: Conformational coupling in ion channels: Conformational coupling in BK potassium channels
.
J. Gen. Physiol.
140
:
625
634
.
Horrigan
F.T.
,
Aldrich
R.W.
.
2002
.
Coupling between voltage sensor activation, Ca2+ binding and channel opening in large conductance (BK) potassium channels
.
J. Gen. Physiol.
120
:
267
305
.
Hoshi
T.
,
Zagotta
W.N.
,
Aldrich
R.W.
.
1994
.
Shaker potassium channel gating. I: Transitions near the open state
.
J. Gen. Physiol.
103
:
249
278
.
Im
W.
,
Roux
B.
.
2001
.
Brownian dynamics simulations of ions channels: A general treatment of electrostatic reaction fields for molecular pores of arbitrary geometry
.
J. Chem. Phys.
115
:
4850
4861
.
Im
W.
,
Roux
B.
.
2002
.
Ion permeation and selectivity of OmpF porin: a theoretical study based on molecular dynamics, Brownian dynamics, and continuum electrodiffusion theory
.
J. Mol. Biol.
322
:
851
869
.
Islas
L.D.
,
Sigworth
F.J.
.
1999
.
Voltage sensitivity and gating charge in Shaker and Shab family potassium channels
.
J. Gen. Physiol.
114
:
723
742
.
Jensen
M.Ø.
,
Jogini
V.
,
Borhani
D.W.
,
Leffler
A.E.
,
Dror
R.O.
,
Shaw
D.E.
.
2012
.
Mechanism of voltage gating in potassium channels
.
Science.
336
:
229
233
.
Karlin
A.
1967
.
On the application of “a plausible model” of allosteric proteins to the receptor for acetylcholine
.
J. Theor. Biol.
16
:
306
320
.
Kohen
D.
,
Tannor
D.
,
.
2000
.
Phase space approach to dissipative molecular dynamics
.
In
Vol. 111
.
Prigogine
I.
,
Rice
S.A.
,
editors
.
John Wiley & Sons, Inc.
,
New York, NY, USA
.
219
398
.
Kramers
H.A.
1940
.
Brownian motion in a field of force and the diffusion model of chemical reactions
.
Physica.
7
:
284
304
.
Kubo
R.
1966
.
The fluctuation-dissipation theorem
.
Rep. Prog. Phys.
29
:
255
284
.
Lange
O.F.
,
Grubmüller
H.
.
2006
.
Collective Langevin dynamics of conformational motions in proteins
.
J. Chem. Phys.
124
:
214903
214923
.
Ledwell
J.L.
,
Aldrich
R.W.
.
1999
.
Mutations in the S4 region isolate the final voltage-dependent cooperative step in potassium channel activation
.
J. Gen. Physiol.
113
:
389
414
.
Leffler
J.E.
1953
.
Parameters for the description of transition states
.
Science.
117
:
340
341
.
Liebovitch
L.S.
1989
.
Analysis of fractal ion channel gating kinetics: kinetic rates, energy levels, and activation energies
.
Math. Biosci.
93
:
97
115
.
Long
S.B.
,
Campbell
E.B.
,
Mackinnon
R.
.
2005
.
Crystal structure of a mammalian voltage-dependent Shaker family K+ channel
.
Science.
309
:
897
903
.
Lopes
P.E.M.
,
Roux
B.
,
Mackerell
A.D.
Jr
.
2009
.
Molecular modeling and dynamics studies with explicit inclusion of electronic polarizability. Theory and applications
.
Theor. Chem. Acc.
124
:
11
28
.
Mannuzzu
L.M.
,
Moronne
M.M.
,
Isacoff
E.Y.
.
1996
.
Direct physical measure of conformational rearrangement underlying potassium channel gating
.
Science.
271
:
213
216
.
Marks
T.N.
,
Jones
S.W.
.
1992
.
Calcium currents in the A7r5 smooth muscle-derived cell line. An allosteric model for calcium channel activation and dihydropyridine agonist action
.
J. Gen. Physiol.
99
:
367
390
.
McManus
O.B.
,
Weiss
D.S.
,
Spivak
C.E.
,
Blatz
A.L.
,
Magleby
K.L.
.
1988
.
Fractal models are inadequate for the kinetics of four different ion channels
.
Biophys. J.
54
:
859
870
.
Mel’nikov
V.I.
,
Meshkov
S.V.
.
1986
.
Theory of activated rate processes: Exact solution of the Kramers problem
.
J. Chem. Phys.
85
:
1018
1027
.
Meyer
R.
,
Heinemann
S.H.
.
1997
.
Temperature and pressure dependence of Shaker K+ channel N- and C-type inactivation
.
Eur. Biophys. J.
26
:
433
445
.
Millhauser
G.L.
,
Salpeter
E.E.
,
Oswald
R.E.
.
1988
.
Diffusion models of ion-channel gating and the origin of power-law distributions from single-channel recording
.
85
:
1503
1507
.
Monod
J.
,
Wyman
J.
,
Changeux
J.P.
.
1965
.
On the nature of allosteric transitions: a plausible model
.
J. Mol. Biol.
12
:
88
118
.
Murata
Y.
,
Iwasaki
H.
,
Sasaki
M.
,
Inaba
K.
,
Okamura
Y.
.
2005
.
Phosphoinositide phosphatase activity coupled to an intrinsic voltage sensor
.
Nature.
435
:
1239
1243
.
Muroi
Y.
,
Arcisio-Miranda
M.
,
Chowdhury
S.
,
Chanda
B.
.
2010
.
Molecular determinants of coupling between the domain III voltage sensor and pore of a sodium channel
.
Nat. Struct. Mol. Biol.
17
:
230
237
.
Neher
E.
,
Sakmann
B.
.
1976
.
Single-channel currents recorded from membrane of denervated frog muscle fibres
.
Nature.
260
:
799
802
.
Neher
E.
,
Stevens
C.F.
.
1977
.
Conductance fluctuations and ionic pores in membranes
.
Annu. Rev. Biophys. Bioeng.
6
:
345
381
.
Niss
M.
2005
.
History of the Lenz-Ising model 1920–1950: from ferromagnetic to cooperative phenomena
.
Arch. Hist. Exact Sci.
59
:
267
318
.
Noceti
F.
,
Baldelli
P.
,
Wei
X.
,
Qin
N.
,
Toro
L.
,
Birnbaumer
L.
,
Stefani
E.
.
1996
.
Effective gating charges per channel in voltage-dependent K+ and Ca2+ channels
.
J. Gen. Physiol.
108
:
143
155
.
Osteen
J.D.
,
Barro-Soria
R.
,
Robey
S.
,
Sampson
K.J.
,
Kass
R.S.
,
H.P.
.
2012
.
Allosteric gating mechanism underlies the flexible gating of KCNQ1 potassium channels
.
109
:
7103
7108
.
Ostmeyer
J.
,
Chakrapani
S.
,
Pan
A.C.
,
Perozo
E.
,
Roux
B.
.
2013
.
Recovery from slow inactivation in K+ channels is controlled by water molecules
.
Nature.
501
:
121
124
.
Pantazis
A.
,
Gudzenko
V.
,
Savalli
N.
,
Sigg
D.
,
Olcese
R.
.
2010
.
Operation of the voltage sensor of a human voltage- and Ca2+-activated K+ channel
.
107
:
4459
4464
.
Papazian
D.M.
,
Shao
X.M.
,
Seoh
S.A.
,
Mock
A.F.
,
Huang
Y.
,
Wainstock
D.H.
.
1995
.
Electrostatic interactions of S4 voltage sensor in Shaker K+ channel
.
Neuron.
14
:
1293
1301
.
Pathak
M.
,
Kurtz
L.
,
Tombola
F.
,
Isacoff
E.
.
2005
.
The cooperative voltage sensor motion that gates a potassium channel
.
J. Gen. Physiol.
125
:
57
69
.
Pathak
M.M.
,
Yarov-Yarovoy
V.
,
Agarwal
G.
,
Roux
B.
,
Barth
P.
,
Kohout
S.
,
Tombola
F.
,
Isacoff
E.Y.
.
2007
.
Closing in on the resting state of the Shaker K+ channel
.
Neuron.
56
:
124
140
.
Perozo
E.
,
Cortes
D.M.
,
Cuello
L.G.
.
1999
.
Structural rearrangements underlying K+-channel activation gating
.
Science.
285
:
73
78
.
Pollak
E.
1986
.
Theory of activated rate processes: A new derivation of Kramers’ expression
.
J. Chem. Phys.
85
:
865
867
.
Qin
F.
2007
.
Principles of single-channel kinetic analysis
.
Methods Mol. Biol.
403
:
253
286
.
Rothberg
B.S.
,
Magleby
K.L.
.
2000
.
Voltage and Ca2+ activation of single large-conductance Ca2+-activated K+ channels described by a two-tiered allosteric gating mechanism
.
J. Gen. Physiol.
116
:
75
100
.
Roux
B.
1995
.
The calculation of the potential of mean force using computer simulations
.
Comput. Phys. Commun.
91
:
275
282
.
Roux
B.
1997
.
Influence of the membrane potential on the free energy of an intrinsic protein
.
Biophys. J.
73
:
2980
2989
.
Roux
B.
1999
.
Statistical mechanical equilibrium theory of selective ion channels
.
Biophys. J.
77
:
139
153
.
Roux
B.
2008
.
The membrane potential and its representation by a constant electric field in computer simulations
.
Biophys. J.
95
:
4205
4216
.
Roux
B.
,
Karplus
M.
.
1991
.
Ion transport in a gramicidin-like channel: dynamics and mobility
.
J. Phys. Chem.
95
:
4856
4868
.
Rowley
C.N.
,
Roux
B.
.
2013
.
A computational study of barium blockades in the KcsA potassium channel based on multi-ion potential of mean force calculations and free energy perturbation
.
J. Gen. Physiol.
142
:
451
463
.
Ryu
S.
,
Yellen
G.
.
2012
.
Charge movement in gating-locked HCN channels reveals weak coupling of voltage sensors and gate
.
J. Gen. Physiol.
140
:
469
479
.
Sasaki
M.
,
Takagi
M.
,
Okamura
Y.
.
2006
.
A voltage sensor-domain protein is a voltage-gated proton channel
.
Science.
312
:
589
592
.
Schnakenberg
J.
1976
.
Network theory of microscopic and macroscopic behavior of master equation systems
.
Rev. Mod. Phys.
48
:
571
585
.
Schoppa
N.E.
,
Sigworth
F.J.
.
1998a
.
Activation of shaker potassium channels. I. Characterization of voltage-dependent transitions
.
J. Gen. Physiol.
111
:
271
294
.
Schoppa
N.E.
,
Sigworth
F.J.
.
1998b
.
Activation of Shaker potassium channels. III. An activation gating model for wild-type and V2 mutant channels
.
J. Gen. Physiol.
111
:
313
342
.
Schreiber
G.
,
Fersht
A.R.
.
1995
.
Energetics of protein-protein interactions: Analysis of the barnase-barstar interface by single mutations and double mutant cycles
.
J. Mol. Biol.
248
:
478
486
.
Schreiber
M.
,
Yuan
A.
,
Salkoff
L.
.
1999
.
Transplantable sites confer calcium sensitivity to BK channels
.
Nat. Neurosci.
2
:
416
421
.
Seoh
S.A.
,
Sigg
D.
,
Papazian
D.M.
,
Bezanilla
F.
.
1996
.
Voltage-sensing residues in the S2 and S4 segments of the Shaker K+ channel
.
Neuron.
16
:
1159
1167
.
T.
,
Yifrach
O.
.
2013
.
Using hierarchical thermodynamic linkage analysis to study ion channel gating
.
J. Gen. Physiol.
141
:
507
510
.
Shrivastava
I.H.
,
Sansom
M.S.
.
2000
.
Simulations of ion permeation through a potassium channel: molecular dynamics of KcsA in a phospholipid bilayer
.
Biophys. J.
78
:
557
570
.
Sigg
D.
2013
.
A linkage analysis toolkit for studying allosteric networks in ion channels
.
J. Gen. Physiol.
141
:
29
60
.
Sigg
D.
,
Bezanilla
F.
.
1997
.
Total charge movement per channel. The relation between gating charge displacement and the voltage sensitivity of activation
.
J. Gen. Physiol.
109
:
27
39
.
Sigg
D.
,
Bezanilla
F.
.
2003
.
A physical model of potassium channel activation: from energy landscape to gating kinetics
.
Biophys. J.
84
:
3703
3716
.
Sigg
D.
,
Stefani
E.
,
Bezanilla
F.
.
1994
.
Gating current noise produced by elementary transitions in Shaker potassium channels
.
Science.
264
:
578
582
.
Sigg
D.
,
Qian
H.
,
Bezanilla
F.
.
1999
.
Kramers’ diffusion theory applied to gating kinetics of voltage-dependent ion channels
.
Biophys. J.
76
:
782
803
.
Sigg
D.
,
Bezanilla
F.
,
Stefani
E.
.
2003
.
Fast gating in the Shaker K+ channel and the energy landscape of activation
.
100
:
7611
7615
.
Sigworth
F.J.
1980
.
The variance of sodium current fluctuations at the node of Ranvier
.
J. Physiol.
307
:
97
129
.
Sigworth
F.J.
1981
.
Covariance of nonstationary sodium current fluctuations at the node of Ranvier
.
Biophys. J.
34
:
111
133
.
Smith-Maxwell
C.J.
,
Ledwell
J.L.
,
Aldrich
R.W.
.
1998
.
Uncharged S4 residues and cooperativity in voltage-dependent potassium channel activation
.
J. Gen. Physiol.
111
:
421
439
.
Soler-Llavina
G.J.
,
Chang
T.-H.
,
Swartz
K.J.
.
2006
.
Functional interactions at the interface between voltage-sensing and pore domains in the Shaker KV channel
.
Neuron.
52
:
623
634
.
Straub
J.E.
,
Borkovec
M.
,
Berne
B.J.
.
1987
.
Calculation of dynamic friction on intramolecular degrees of freedom
.
J. Phys. Chem.
91
:
4995
4998
.
Tiwari-Woodruff
S.K.
,
Schulteis
C.T.
,
Mock
A.F.
,
Papazian
D.M.
.
1997
.
Electrostatic interactions between transmembrane segments mediate folding of Shaker K+ channel subunits
.
Biophys. J.
72
:
1489
1500
.
Tiwari-Woodruff
S.K.
,
Lin
M.A.
,
Schulteis
C.T.
,
Papazian
D.M.
.
2000
.
Voltage-dependent structural interactions in the Shaker K+ channel
.
J. Gen. Physiol.
115
:
123
138
.
Tolokh
I.S.
,
White
G.W.N.
,
Goldman
S.
,
Gray
C.G.
.
2002
.
Prediction of ion channel transport from Grote–Hynes and Kramers theories
.
Mol. Phys.
100
:
2351
2359
.
Torrie
G.M.
,
Valleau
J.P.
.
1977
.
Nonphysical sampling distributions in Monte Carlo free-energy estimation: Umbrella sampling
.
J. Comput. Phys.
23
:
187
199
.
Van Kampen
N.G.
1992
.
Stochastic processes in physics and chemistry
. 2nd ed.
North-Holland
,
Amsterdam
.
480 pp
.
Vandenberg
C.A.
,
Bezanilla
F.
.
1991
.
A sodium channel gating model based on single channel, macroscopic ionic, and gating currents in the squid giant axon
.
Biophys. J.
60
:
1511
1533
.
Vargas
E.
,
Yarov-Yarovoy
V.
,
Khalili-Araghi
F.
,
Catterall
W.A.
,
Klein
M.L.
,
Tarek
M.
,
Lindahl
E.
,
Schulten
K.
,
Perozo
E.
,
Bezanilla
F.
,
Roux
B.
.
2012
.
An emerging consensus on voltage-dependent gating from computational modeling and molecular dynamics simulations
.
J. Gen. Physiol.
140
:
587
594
.
Waldram
J.R.
1985
.
The Theory of Thermodynamics
.
Cambridge University Press
,
Cambridge
.
352 pp
.
White
G.W.N.
,
Goldman
S.
,
Gray
C.G.
.
2000
.
Test of rate theory transmission coefficient algorithms. An application to ion channels
.
Mol. Phys.
98
:
1871
1885
.
Wyman
J.
1967
.
.
J. Am. Chem. Soc.
86
:
2202
2218
.
Yang
N.
,
Horn
R.
.
1995
.
Evidence for voltage-dependent S4 movement in sodium channels
.
Neuron.
15
:
213
218
.
Yarov-Yarovoy
V.
,
DeCaen
P.G.
,
Westenbroek
R.E.
,
Pan
C.-Y.
,
Scheuer
T.
,
Baker
D.
,
Catterall
W.A.
.
2012
.
Structural basis for gating charge movement in the voltage sensor of a sodium channel
.
109
:
E93
E102
.
Yifrach
O.
,
MacKinnon
R.
.
2002
.
Energetics of pore opening in a voltage-gated K+ channel
.
Cell.
111
:
231
239
.
Zagotta
W.N.
,
Hoshi
T.
,
Aldrich
R.W.
.
1994
.
Shaker potassium channel gating. III: Evaluation of kinetic models for activation
.
J. Gen. Physiol.
103
:
321
362
.
Zandany
N.
,
M.
,
Orr
I.
,
Yifrach
O.
.
2008
.
Direct analysis of cooperativity in multisubunit allosteric proteins
.
105
:
11697
11702
.

Abbreviations used in this paper:

• DSM

discrete-state Markov

•
• GH

Grote–Hynes

•
• GLE

generalized Langevin equation

•
• HH

Hodgkin–Huxley

•
• LF

large friction

•
• LFER

linear free energy relationship

•
• MD

molecular dynamics

•
• MWC

Monod–Wyman–Changeux

•
• PMF

potential of mean force

•
• TST

transition state theory

•
• VSD

voltage-sensing domain