INTRODUCTION

TOP ABSTRACT INTRODUCTION MATERIALS AND METHODS RESULTS AND DISCUSSION CONCLUSIONS REFERENCES

Determination of the crystal structure of the KcsA potassium channel (Doyle et al., 1998) and its subsequent refinement at 2.0-Å resolution (Morais-Cabral et al., 2001; Zhou et al., 2001) have provided valuable information on this biologically important class of channels, especially on the operation of the selectivity filter. There are many different types of potassium channels, which differ widely in their conductances and gating characteristics (Hille, 2001). Some of these are voltage gated or activated by intracellular Ca²⁺ ions, whereas the activation of the bacterial KcsA channel depends sensitively on internal pH (Cuello et al., 1998; Heginbotham et al., 1999). Conductance levels of various types of potassium channels range from 4 to 270 pS (for example, see Dittrich and Daut, 1999; Dopico et al., 1999; Hirsch et al., 1999; Hirschberg et al., 1999; Lara et al., 1999; Latorre et al., 1989; Noulin et al., 1999; Reid et al., 1999). Despite this diversity, they all share the common feature of being highly selective to potassium ions (Hille, 2001) and display broadly similar selectivity sequences and relative permeability ratios for monovalent cations (Christophersen, 1991; Hill et al., 1989; Hu et al., 1989; Gelband and McCullough, 1993; LeMasurier et al., 2001; Schlatter et al., 1993; Shapiro and DeCoursey, 1991; Shvinka and Caffier, 1988; Tabcharani and Misler, 1989).

The TVGYGD amino acid sequence of the peptide chains lining the selectivity filter of all potassium channels is known to be highly conserved (Heginbotham et al., 1992; Schrempf et al., 1995; Doyle et al., 1998, MacKinnon et al., 1998). The only charged residue in the external pore region of KcsA that is not conserved is the glutamic acid E₇₁. Because this residue appears to be protonated (Ranatunga et al., 2001), its presence is not essential to channel function. Thus, it is likely that the diversity of potassium channels results from structural changes on the protein architecture near the intracellular segment of the pore, which have very different sequences. While the search for the complete tertiary structure of potassium channels continues (for example, see Hong and Miller 2000; Li-Smerin et al., 2000; Perozo, 2000; Lu et al., 2001), useful insights to the structure of the pore may be obtained from a study of the inverse problem, that is, predicting relevant aspects of the channel structure from its functional properties.

Appearance of the tertiary structure of the KcsA protein has stimulated much activity in the modeling of ion permeation in potassium channels (for reviews, see Roux et al., 2000; Kuyucak et al., 2001; Tieleman et al., 2001). Most of these studies focus on the selectivity filter and attempt to understand its operational principles by performing molecular dynamics (MD) simulations (Åqvist and Luzhkov, 2000; Guidoni et al., 2000; Shrivastava and Sansom, 2000; Allen et al., 2000; Bernéche and Roux, 2001). However, MD simulations are too slow at present to determine the channel conductance, which is the most important functional property of an ion channel. To facilitate the calculation of conductance, one has to use a more coarse-grained simulation method such as Brownian dynamics (BD) or Poisson-Nernst-Planck equations. In view of the difficulties associated with the application of continuum theories of electro-diffusion to narrow pores (Corry et al., 2000; Graf et al., 2000), BD appears to be the method of choice for this purpose (Bek and Jakobsson, 1994; Chung et al., 1999, 2002; Allen and Chung 2001; Mashl et al., 2001).

In our previous studies of the KcsA channel (Allen and Chung, 2001; Chung et al., 1999, 2002), we investigated its permeation properties by calculating the electrostatic potential energy profiles for multiple K⁺ ions and performing three-dimensional BD simulations. Electrostatic calculations show that, without the stabilizing effect of backbone dipole and charged side-chains in the protein, a potassium ion attempting to cross-the channel would face an insurmountable energy barrier. When these charges are placed inside the channel wall, this large barrier turns into a deep potential well, which accommodates two K⁺ ions in the selectivity filter and one in the cavity in the absence of an applied potential. By simulating the trajectories of ions using BD, we deduced many of the experimentally observed properties of the channel. Among these are the channel conductance, the current-voltage relationships, the conductance-concentration curves, and the reversal potentials with asymmetrical ionic concentrations in the two sides of the channel.

Here we use BD simulations to explore whether the widely differing properties of potassium channels found in nature can be understood by small modifications of the channel geometry. For this purpose, we first construct a conducting state of the KcsA potassium channel that includes all the experimentally determined channel protein. From this detailed model, we extract a simplified KcsA model by replacing the atoms in the protein with a homogeneous, low dielectric medium, and representing the polar groups and oppositely charged pairs, deemed responsible for ion permeation, as dipoles. Also the irregular water-protein interface found in the crystal structure is smoothed out. The atomic-detail and simplified models of KcsA are compared to ascertain that they have qualitatively similar functional properties. These simplifications reduce the computational time for BD simulations, which is important for a time consuming systematic study. But more importantly, the simplified structure provides a template for modeling of the diverse range of potassium channels using BD. Using a simplified structure for this purpose (rather than the atomic structure of KcsA) is more sensible because, as noted above, there are large variations in their sequences. Here, we explore the influence of the intrapore geometry on the diversity question by changing the radius of the pore entrance on the intracellular side systematically from 1.5 to 5 Å while keeping the dimensions of the selectivity filter and cavity constant. We make several predictions about the properties of the potassium channels with large and small conductances.

It is worthwhile to emphasize that such a procedure would be much harder to use with the MD method, even if the current time constraint for the calculation of conductance were not a problem. This is because MD is much more sensitive to details in channel structure, and finding the correct atomic structure that would reproduce the properties of the channel would be a very time consuming task. We emphasize that this is simply a phase-space problem in the inverse method (i.e., going from function to structure) and has nothing to do with the capabilities of MD. Thus MD studies of eukaryotic potassium channels have to wait for the solution of their crystal structures, which may take a long time. In contrast, the implicit treatment of water and averaging over the atomic structure of protein in BD renders it insensitive to such details, making it an ideal tool for exploring the structure-function relationships in a diverse range of potassium channels.

	MATERIALS AND METHODS

TOP ABSTRACT INTRODUCTION MATERIALS AND METHODS RESULTS AND DISCUSSION CONCLUSIONS REFERENCES

Channel models

Electron paramagnetic resonance studies with site-directed spin labeling (Perozo et al., 1998, 1999; Gross et al., 1999; Liu et al., 2001) have indicated that the crystal structure of KcsA (Doyle et al., 1998) corresponds to the closed state of the channel, and the transmembrane helices lining the intracellular pore (TM2) move away from the channel axis during gating. To study the conductance properties of KcsA, therefore, we need to construct an open-state model by increasing the intrapore radius of the crystal structure (PDB structure 1BL8). This is achieved by carrying out MD simulations with the CHARMM package using the CHARMM 19 extended atom parameter set, as described previously (Allen et al., 2000). A reservoir containing extended simple point charge water molecules is attached to each end of the experimentally determined protein, and the pore formed by the subunits is allowed to be hydrated slowly for 200 ps. Then, three K⁺ ions are placed near the positions determined by Rb⁺ difference maps (Doyle et al., 1998) and held with 100 kT/Ų harmonic constraints. An additional K⁺ ion is placed near the intracellular entrance of the channel. Then, in a series of MD simulations, each lasting 100 ps, atomic constraints of 20 kT/Ų are applied to -carbon atoms that moves them by 0.1 Å from their current positions, and these target restraints are moved away from the channel axis until the minimal radius of the intracellular pore is ~3 Å. Finally, all the restraints on the inner and outer helices are removed so that transmembrane helices are free to rotate during a further 50 ps of simulation. Water molecules and ions are discarded, and a single subunit is chosen and replicated to impose a fourfold symmetry about the channel axis.

Detailed channel model

View larger version (45K): [in this window] [in a new window]   FIGURE 1   Models of potassium channels. (A) Shape of the KcsA channel is modified such that the minimal radius of the wider segment of the the pore is 3 Å. The dark gray, white, black, and light gray ribbons represent the outer helices, pore helices, selectivity filter, and inner helices, respectively. (B) Solid line shows the outline of a simplified model channel with the radius of the wider segment. The shape of this channel closely corresponds to that of the modified KcsA channel shown in A. The positions of dipoles on the channel wall are indicated: , 10 of the 20 carbonyl and hydroxyl oxygen atoms; , N termini of the helix dipoles; and , mouth dipoles. (C) Outlines of a set of model channels are superimposed. The radius of the channel facing the intracellular space (left side) is systematically changed from 1.5 to 5 Å.

Simplified channel model

30

Reservoirs

Electrostatic calculations

The selectivity filter and cavity of the potassium channel are usually occupied by two to three K⁺ ions. When an ion in the cavity attempts to move toward the intracellular or extracellular side, it encounters an energy barrier. This barrier plays a pivotal role in the permeation process, and therefore, we wish to visualize it by constructing potential energy profiles in the presence of several K⁺ ions in the channel with different strengths of the applied electric field. We construct these profiles by moving a test ion in 1-Å intervals, holding it fixed at each step. We then allow the other ions resident in the channel to adjust their positions so that the force on them vanishes, thus minimizing the total energy of the system. The minimization is performed at each step, and the positions of the resident ions and the potential energy of the test ion are recorded. This energy corresponds to the work required to bring in the test ion from an infinite distance to its current position, and it is given by

(1)

in which the index i refers to the test ion and j ranges over all the other ions. The first four terms in Eq. 1 represent, respectively: 1) the external potential due to the applied field, fixed charges, and charges induced by these; 2) the self potential due to the surface charges induced by the test ion on the channel boundary; 3) the image potential due to the charges induced by the ion j; and 4) the Coulomb potential due to the ion j. The last term, U_B,i, represents the Born energy for an ion entering from the reservoir to the channel, whose dielectric constants are assumed to be  = 60 for water and 2 for the protein. The reasons for choosing these values of  are given in Chung et al. (1999, 2002). This Born energy barrier of height of 0.6 kT is implemented as described in Chung et al. (1999).

Brownian dynamics

We use BD simulations to determine the channel conductance under various conditions. The algorithm for BD is conceptually simple: the motion of an ion with mass m_i and charge q_i in a fluid is governed by the Langevin equation

(2)

By solving the Langevin equation at discrete time steps for each ion in the system, we follow the trajectories of all the ions in the assembly and determine the number of ions crossing the channel in a given simulation period. The four terms on the right-hand side of Eq. 2 correspond to, respectively, an average frictional force with a friction coefficient m_ii, the stochastic force arising from random collisions with water molecules, the electric and the short-range forces experienced by the ion (Chung et al., 2002). The electric field, E, is obtained from the numerical solution of Poisson's equation on a grid of points and stored in a number of lookup tables (Hoyles et al., 1998b). During simulations, the field and potential at desired positions are reconstructed by interpolating between the table entries. The method, described and validated in Hoyles et al. (1998a), drastically reduces the computational time and thus enabling computation of the conductance of model channels. For further technical details and the computational steps involved in the implementation of the BD algorithm, we refer to Chung et al. (1998, 1999) and Corry et al. (2001), where the method is applied to model nicotinic acetylcholine receptor, potassium, and calcium channels, respectively.

We generate the potential difference across the channel via an applied electric field of constant strength E, which could be due to a pair of isopotential voltage plates located far from the channel. Assuming that the space outside the reservoirs is composed of an electrolyte solution, the potential drop due to this field occurs mainly across the membrane and the channel. In the absence of any dielectric boundary, the potential difference across a channel of length l would be simply lE. The presence of a dielectric boundary and charges on the protein wall, however, severely distorts the applied field. Thus, the potential difference across the channel depends on the selected reference points at the two sides of the channel. A convenient choice for this purpose is the center of each reservoir at z = ±50 Å, which we use throughout this work. Due to the asymmetric distribution of charge residues, there is a potential difference of V = 23 mV between these two points even in the absence of an applied field. When a uniform electric field of strength E = n × 10⁷ V/m is applied, the potential difference between the reference points becomes V = 23 + n × 134 mV. Thus, for the range of applied fields considered in this work with n = 2, 1, 1, 2, the corresponding potentials are V = 291, 157, 111, 245 mV, respectively.

Bulk ionic diffusion coefficients are used everywhere except in the selectivity filter, where it is reduced to 50% of the bulk value for K⁺ ions. As shown in Chung et al. (1999), BD results are insensitive to the precise values of the diffusion coefficients, especially in the filter region. With the exception of the current-concentration curves, we carry out the simulations with a total of 64 ions in the system (32 K⁺ and 32 Cl), corresponding to the physiological concentration of 300 mM KCl when the length of each reservoir is fixed at 32.5 Å. This higher concentration is used so as to obtain better statistics. We consider only symmetric solutions in this work; hence the ions are distributed evenly between the intracellular and extracellular reservoirs.

The BD program is written in FORTRAN following the algorithm of van Gunsteren and Berendsen (1982), vectorized and executed on a supercomputer (Fujitsu VPP-300). The amount of vectorization varies from 67% to 92%, depending on the number of ions in the reservoirs. With 64 ions, the CPU time needed to complete a simulation period of 1 µs (10 million time steps) is ~28 h for simplified channels and 53 h for the channel with the full atomic details.

Throughout we quote energy in room temperature units (kT) and dipole moment in Coulomb-meter (Cm). We note 1 kT = 4.11 × 10²¹ J or 0.592 Kcal/mol and 1 Debye = 3.33 × 10³⁰ Cm. The following physical constants are used in our calculations: (Note that the friction coefficient is related to the diffusion coefficient via the Einstein relation D = kT/m.)

Masses: m_K = 6.5 × 10²⁶ kg,     m_Cl = 5.9 × 10²⁶ kg.

Diffusion coefficients:  D_K = 1.96 × 10⁹ m²/s,
D_Cl = 2.03 × 10⁹ m²/s.

Ionic radii: R_K = 1.33 Å,     R_Cl = 1.81 Å.

Room temperature: T = 298 K.

	RESULTS AND DISCUSSION

TOP ABSTRACT INTRODUCTION MATERIALS AND METHODS RESULTS AND DISCUSSION CONCLUSIONS REFERENCES

Optimization of dipole strengths in detailed model

Of the nine ionizable residues on each subunit of the KcsA protein, we have determined that only the aspartate-arginine pair near the extracellular entrance (D₈₀-R₆₄), and the glutamate-arginine pair near the intracellular entrance (E₁₁₈-R₁₁₇) require charge to enable conduction (Allen and Chung, 2001; Chung et al., 2002). The current across the channel is dependent on the magnitude of charges placed on these four residues.

We assign D₈₀ and R₆₄ charges q_e and +q_e to create four extracellular mouth dipoles, and place charges q_i and +q_i on E₁₁₈ and R₁₁₇ to create four intracellular mouth dipoles. In Fig. 2 A, we show the variation in outward (filled circles) and inward (open circles) currents with charge on the intracellular dipole q_i. These results are obtained from the model with the intrapore radius of 4 Å with a symmetric 300 mM KCl solution and an applied field of ±2 × 10⁷ V/m. Each point in this and subsequent figures, unless otherwise stated, is the average of 10 simulations, each lasting 0.1 µs (or 1 million time steps). The error bar accompanying each data point is one standard error of means and is not shown if it is smaller than the size of the data point. The full charge of ±e is placed on the aspartate-arginine pair during optimization of q_i. Both outward and inward currents increase with q_i reaching a maximum at 0.5 e, and then decrease as the charge is further increased. We next fix q_i = 0.5 e and determine the variation in outward and inward currents with charge on the extracellular dipoles. As shown in Fig. 2 B, both outward (filled circles) and inward (open circles) currents increase steadily with q_e, reaching the maximum when the full unit charge of ±e is placed on the aspartate-arginine pair. From these series of simulations, we assume that the D₈₀ and R₆₄ residues are fully charged in a conducting state (q_e = e), whereas the E₁₁₈ and R₁₁₇ residues guarding the intracellular entrance of the channel are partially charged (q_i = 0.4 e). For the simplified channel models, we place the same charges of q_i = 0.4 e and q_e = e on the intracellular and extracellular dipoles, respectively.

View larger version (13K): [in this window] [in a new window]   FIGURE 2   Changes in channel currents with the strength of mouth dipoles. The current across the channel with full atomic details is plotted against the strength of mouth dipoles lining the intracellular entrance of the channel. Each point in the figure is derived from a simulation period of 1 µs, or ten-million time steps. The error bar in this and all subsequent figures represent 1 SE, and is not shown when it is smaller than the plotting symbols. A uniform electric field of 2 × 107 V/m is applied in the z direction. (A) Outward () and inward () currents are plotted against the charges on the E118 and R117 residues. The full unit charges are placed on the D80 and R64 residues throughout. (B) Outward () and inward () currents are plotted against the charges on the D80 and R64 residues. The charges placed on the E118 and R117 residues are 0.5 e throughout.

Detailed model versus simplified model

Here we show that many of the channel properties obtained from a simplified model are qualitatively similar to those obtained from the model with the full atomic details. These include the potential energy profiles, distributions of ions in the channel, channel conductance, and conductance-concentration curves.

In Fig. 3, we compare the potential energy profiles obtained from the atomic-detail model to those obtained from the simplified model with an intrapore radius of 3.5 Å. (Note that the atomic-detail model has a minimal radius of 3 Å and its average radius is slightly larger than that.) With no ions in the channels, an ion traversing along the central axis encounters a deep energy well, created by the dipoles and ionizable residues. The depth of the energy well reaches to 57 kT for the simplified model (solid line in Fig. 3 A) and 66 kT for the atomic-detail model (broken line). The deep energy well attracts two ions into the selectivity filter, which remain in a stable equilibrium even in the presence of an applied field of 2 × 10⁷ V/m. The presence of the two ions in the selectivity filter drastically alters the profile seen by a third ion approaching the cavity from inside (Fig. 3 B). In both channel models, the third ion sees an energy well with its minimum located at z = 20 Å. From there, the ion has to climb over a shallow barrier of height of 3.2 kT for the detailed model (broken line) and 2.8 kT for the simplified model (solid line) to proceed toward the cavity. Although the wells are deeper in the detailed model compared with the simplified one, their general shapes and the heights of the residual barriers are similar.

View larger version (13K): [in this window] [in a new window]   FIGURE 3   The potential energy profiles obtained from the detailed model (broken line) and the simplified model with the intrapore radius of 3.5 Å (solid lines). The profiles illustrated in A are constructed with no resident ion in the channel, whereas those illustrated in B are constructed with two ions in the selectivity filter, one at z = 15 and the other at z = 21 Å.

To determine the locations of ions, we divide the model channel into 100 thin sections and tabulate the average number of ions over a 0.1-µs simulation. These are illustrated with histograms in Fig. 4 for the simplified model (A) and the detailed model (B). In the presence of an applied field of 2 × 10⁷ V/m, there are two prominent peaks in the selectivity filter and another peak near the intracellular entrance of the pore. On average, there are 2.3 ions in the selectivity filter and cavity, and 0.5 ions (Fig. 4 A) and 1.2 (Fig. 4 B) near the intracellular entrance of the channel. The centers of the maxima for the histograms illustrated in Fig. 4 are at z = 10.7 (near the T₇₅ carbonyl oxygen) and 18.9 (near the Y₇₈ carbonyl oxygen) and z = 21 Å (near the D₈₀ side-chain. The peaks in the selectivity filter are separated by a mean distance of 7.3 Å. We note that, in the absence of an applied field (not shown), there appears an additional peak in the cavity in addition to the two peaks in the selectivity filter, consistent with the observations from the x-ray diffraction maps (Doyle et al., 1998; Morais-Cabral et al., 2001).

View larger version (21K): [in this window] [in a new window]   FIGURE 4   Ions in the channel. The histograms indicate the locations of ions in the simplified model (A) and the detailed model (B). Ions drift across from the intracellular space (left) to the extracellular space (right) under the influence of an applied electric field.

As shown in Fig. 5, the current-concentration curves obtained from the two models have the same shape as those observed experimentally in potassium channels (Hille, 2001). The two curves in Fig. 5 are the outward currents obtained from the detailed model (open circles) and the simplified model with the intrapore radius of 3.5 Å (filled circles). Although the current at any concentration is slightly higher for the detailed channel than the simplified channel, the shape of the curve is similar. The current increases rapidly with an increasing ionic concentration initially and then increases at a slower rate with a further increase in concentration. The solid lines drawn through the data points are calculated from the Michaelis-Menten form (Hille 2001).

View larger version (13K): [in this window] [in a new window]   FIGURE 5   Current-concentration relationships. The outward currents across the simplified model () and the detailed model () are plotted against concentrations. The points are fitted by the solid lines using the Michaelis-Menten equation. Each point is derived from the simulation period of 3.2 µs () or 1.3 µs ().

From these comparisons, we conclude that some of the distinctive properties deduced from the channel model incorporating the full atomic details of the protein are broadly similar to those obtained from a simplified model, where the atoms forming the protein are represented as a homogeneous dielectric medium and the ionizable residues as dipoles. Despite the approximations imposed on the simplified model, BD captures the salient conduction properties of potassium channels. This is because the essential features that govern the permeation of ions across a narrow pore are captured in our simplified model. These are the overall channel shape, the magnitude of charges that create an energy well deep enough to attract three ions in the equilibrium condition and two ions in the conducting state, the radius of the intracellular pore and the charges placed on the ionizable residues guarding the channel entrances. In the following simulations, we make use of several simplified channel models to determine how conductance changes with the intrapore radius.

Radius dependence of conductance

Based on intuition gathered from a simple application of Ohm's law or from more sophisticated continuum theories of electrodiffusion, one would expect that the channel conductance is essentially determined by the narrow selectivity filter region, and variations in the pore size near the intracellular mouth would cause relatively small changes in the overall conductance of the channel. To test this notion, we measure currents across seven channel models with different radii of the intrapore entrances. As illustrated in Fig. 6, both outward (Fig. 6 A) and inward (Fig. 6 B) currents across the channel increase steeply with increasing intrapore radius of the channel. With an applied electric field of 2 × 10⁷ V/m and an ionic concentration of 300 mM in the reservoirs, the outward current increases from 0.18 to 48 pA as the radius increases from 2 to 5 Å. When the direction of the field is reversed, the inward current increases from 2.1 to 65 pA as the radius changes from 2.5 to 5 Å (no inward current could be recorded when the intrapore radius is 2 Å.) With an applied field of 10⁷ V/m, the outward and inward currents with the 5-Å channel are 21 pA (187 pS) and 33 pA (210 pS), respectively. At all applied potentials, the channel current is seen to be very sensitive to the intrapore radius, increasing sharply as the intracellular mouth is made wider.

View larger version (9K): [in this window] [in a new window]   FIGURE 6   Dependence of channel currents on the intrapore radius of the channel. The outward (A) and inward (B) currents are plotted against the radius of the intracellular aspect of channel entrance. The applied field to obtain the current is ±2 × 107 V/m.

To understand how this feature arises from the ion-ion and ion-channel interactions, we study in the following potential energy profiles for multiple ions and distribution of ion positions obtained from an analysis of BD simulations. Because of the asymmetric nature of the potassium channel, the outward and inward currents will be considered separately.

Effect of radius on outward currents

In the absence of other ions, there is a very deep potential well in the selectivity filter and no conduction is possible (Chung et al., 1999). Even when there are two ions in the channel, conduction rarely occurs, and the presence of a third ion is required for an optimal operation of the channel. Thus, to elucidate the relationship between the current and the intrapore radius, we examine the potential energy profile of an ion given that the channel is already occupied by two other ions. These potential energy profiles, obtained in the presence of an electric field of 2 × 10⁷ V/m, are illustrated in Fig. 7 A for three channels with radii 2, 3, and 4 Å. Two prominent features of these profiles are the energy barrier centered at z = 10 Å and the accompanying well at z = 20 Å. The height of the energy barrier U an ion needs to surmount to traverse the channel from left to right decreases progressively from 7.7 to 0.8 kT as the radius of the intrapore gate is widened from 2 to 5 Å (Fig. 7 B). The reduction in the barrier height is expected to make an ion's permeation from the inner well to the cavity easier. In the inset, the log of the current across the channel is plotted against U, which clearly shows the exponential decrease in conductance with the increasing barrier height. The other factor that contributes to the increase in conductivity with radius is the probability of an ion entering the channel, which is roughly proportional to the cross-sectional area of the channel entrance.

View larger version (14K): [in this window] [in a new window]   FIGURE 7   Rate-limiting factors in outward currents. (A) Potential energy profiles encountered by an ion traversing along the central axis of the channel when there are two other ions in or prenear the selectivity filter are shown for the channels with radii 2 Å (solid line), 3 Å (long-dashed line), and 4 Å (dashed line). Ions need to climb over the energy barrier, whose height is denoted as U, to move across the channel. (B) The barrier height U is plotted against the intrapore radius. In the inset, the outward current in the logarithmic scale is plotted against U.

The electrostatic and geometric pictures presented above suggest possible mechanisms for the steep increase in channel current with intrapore radius but they do not give the full picture underlying the dynamics of ion permeation, and in particular, explain why the selectivity filter does not create a bottle neck in ion diffusion. For that purpose, one needs a quantitative analysis of BD simulations that are quite laborious. An alternative approach is to create animated visualizations of the BD simulations, which provide immediate insights into the ion permeation problem in potassium channels. Animations clearly show that the trigger of a conduction event is the passage of an ion over the potential barrier at z = 10 Å when there are two ions resident in the selectivity filter. Once an ion is over this barrier, it moves toward the selectivity filter at a speed that is many times greater than the drift velocity in a bulk solution. The reason for this fast diffusion of the third ion in the region between the barrier and the selectivity filter is the enhancement of the applied electric field by the local charges (up to two orders of magnitude for a single ion). The deep energy well created by the carbonyl groups and the macrodipoles are virtually eliminated by the presence of two resident ions. As the third ion approaches the selectivity filter, it invokes a large Coulomb repulsion on the outermost ion, which is expelled from the selectivity filter almost simultaneously. While MD simulations indicate a reduced diffusion coefficient of K⁺ ions in the selectivity filter (Allen et al., 2000), this does not hinder the motion of the third ion much, because it is already near the pore mouth and needs to move only a few angstroms to complete the conduction event. This billiard ball-like mechanism, together with the enhancement of local forces on ions, explain why the selectivity filter does not form a bottle neck. Inferences made from the analysis of the most recent experiments on ion conduction across the selectivity filter of KcsA (Morais-Cabral et al., 2001) strongly support the above permeation mechanism.

This description of ion conduction in the potassium channels can be made more quantitative by analyzing the average times an ion spends in various parts of the channel. The time the energy well remains vacant is clearly dependent on the cross-sectional area of the channel entrance. The time an ion lingers in the well centered at z = 20 Å is expected to be correlated with the barrier height it encounters. We take a snapshot of the channel once every 100 fs for 1 µs and count how many times an ion is present in the intracellular segment between z = 23 and z = 10 Å. When a driving field of 2 × 10⁷ V/m is applied, the average time per conduction event the well remains vacant are 25.5, 4.1, 6.4, 3.4, 1.6, and 1.9 ns for the channels with the intrapore radii of 2.5, 3, 3.5, 4, 4.5, and 5 Å. The corresponding average durations per conduction event that the well remains occupied are 47.2, 8.4, 3.1, 2.9, 2.2, and 1.3 ns. Thus, the time the well remains vacant decreases progressively as the radius is increased, because the probability of an ion entering is approximately proportional to the cross-section of the channel aperture. Once an ion enters the well, the time it remains in the well depends on the height of the barrier, which, as we have shown, decreases with an increasing intrapore radius. These two underlying mechanisms cause a large increase in the channel conduction with a small change in the channel geometry.

Increasing the applied potential decreases both waiting times. Thus, the rate-limiting step for conduction in a channel is the time it takes for an ion to stumble into the energy well created by the mouth dipoles and the time it takes for this ion to climb out of the well. The time it takes for an ion to enter the channel will be reduced if the depth of the energy well is increased (by increasing the dipole strength), but the height of the barrier the ion needs to surmount after it enters the well also increases. For the channel to transfer the maximal number of ions per unit time, the strength of dipoles lining the mouth of the channel must be large to reduce the vacant time but, at the same time, not so large that it renders the potential barrier encountered by an ion insurmountable. This explains why the current in Fig. 2 A decreases after reaching a maximum at q_i = 0.5 e.

To illustrate where in the channel ions reside predominantly, we bisect the channel and tabulate the number of ions in the intracellular or left-hand side segment n_l and in the extracellular or right-hand side segment n_r. The resulting ion configuration is denoted as the state [n_l, n_r]. The most common state of the channel, with an applied field of 2 × 10⁷ V/m is [1, 2], or one ion is present in the intracellular half of the segment and two ions reside in the extracellular half of the segment. The percentage of time the 2-, 3-, 4-, and 5-Å channels are in different states are listed in Table 2. The probability that the channel is in the [1, 2] or [1, 3] state decreases as the intrapore radius increases, whereas the probability that the channel is in the [0, 2] or [0, 3] state increases. Although there are many other transient states, these four states account for more than 90% of the occupation probability.

Effect of radius on inward currents

Although the radius of the selectivity filter and the dipoles lining it and pointing to the cavity remain unchanged, the magnitude of inward currents is also strongly influenced by the intrapore radius (see Fig. 6 B). To explain this seemingly paradoxical situation, we again refer to the multiion potential energy profiles. In Fig. 8 A, we show the energy profiles encountered by a third ion as it moves from the selectivity filter toward the intracellular space (right to left) for the 2.5, 3.5, and 4.5 Å channels. The procedure used to derive these profiles is the same as those used for Fig. 7 A, in that two ions are placed in the selectivity filter and a third ion (the test ion) is moved from z = 10 to z = 40 Å in 1-Å steps. For the test ion to traverse the channel toward the intracellular space, first it has to surmount a large barrier U peaking near z = 10 Å. The ion then lingers on near the intracellular entrance of the channel, where the well created by the E₁₁₈ residues is deepest. Subsequently, the ion is ejected from this well when another ion approaches it. The height of the first barrier decreases with increasing intrapore radius, as shown in Fig. 8 B. In the inset, the log of the inward current is plotted against the height of the barrier U. Fits with a Boltzmann factor indicate a very good correlation between the barrier height U and the channel current. Thus, this barrier seems to be directly responsible for the steep rise in current with the intrapore radius. The increase in the barrier height with decreasing radius is due to the shape of the hydrophobic channel segment changing from a cylinder of a uniform radius to a funnel. The magnitude of induced surface charges on the channel wall, hence the energy barrier, progressively decreases with an increasing aperture of the intracellular side of the channel.

View larger version (15K): [in this window] [in a new window]   FIGURE 8   The rate-limiting factors in inward currents. (A) Potential energy profiles encountered by the test ion traversing along the central axis of the channels with radii 2.5 Å (solid line), 3.5 Å (broken line), and 4.5 Å (dashed line). Two ions are placed in the selectivity filter at z = 14 and z = 18 Å, and a third ion is moved from z = 10 Å toward the intracellular space. (B) Barrier height U is plotted against the intrapore radius. In the inset, the inward current in the logarithmic scale is plotted against U.

The above analysis suggests that the rate-limiting step for conduction is the time it takes for one of the ions located in the selectivity filter/cavity region to break away and climb over the first barrier. Animated visualizations of the BD simulations broadly support this picture. When only two ions are in the selectivity filter, the inner ion can rarely, if ever, break away and traverse the channel. The presence of a third ion in the filter, therefore, is a prerequisite for a conduction event. An analysis of the occupation probabilities of various multiion states given in Table 3 also confirms the above interpretation of the inward conduction events. The channel is seen to be predominantly in the [1, 3] and [0, 3] state, in contrast to the results in Table 2.

View this table: [in this window] [in a new window]   TABLE 3   Occupation probabilities of multiion states [n1, nr] for different channel radii when a field of 2 × 107 V/m is applied in the z direction

Ions in the channel

The specific ion binding sites observed in the x-ray diffraction picture are obtained from a closed channel at liquid nitrogen temperatures (Morais-Cabral et al., 2001). At room temperature, one expects these individual peaks to be washed out and replaced by a broad distribution of ions in the filter (for example, see Fig. 9 A in Chung et al., 2002). Ion distribution in the filter will be further modified in an open, conducting channel. Thus permeation models, especially a coarse-grained model such as BD, cannot be expected to reproduce these binding sites in all microscopic detail. Nevertheless, it is of interest to study the ion distribution in the channel and see whether the average description obtained from the BD simulations is in broad agreement with the experiments. As an illustration of this point, we show in Fig. 9 the distribution of ions within the 3.5-Å channel under the applied field of ±2 × 10⁷ V/m. When the field is in the +z direction pushing ions outward, the resident ions in the selectivity filter tend to dwell around two cites centered at z = 11 and z = 19.3 Å. The average numbers of ions on the right and left half of the channel are 2.1 and 0.5, respectively. The peaks in the histogram represent the locations of the energy minima, which depend on the distribution of the dipoles on the protein wall as well as the direction of the applied electric field. When the direction of the field is reversed so that ions move inward, an additional peak appears near the selectivity filter and inner chamber in the histogram (Fig. 9 B). The peaks in the ion distribution remain at similar positions in the selectivity filter but the ion that dwells preferentially near the inner segment of the selectivity filter now darts back and forth from this position to the inner chamber. The broad, third peak in the histogram is located at z = 6.2 Å. The average number of ions on the right and left half of the channel are 2.3 and 0.3 ions. The small changes in these numbers from the outward conduction case are consistent with the ions waiting on the right side for permeation instead of the left side. The average number of ions in the selectivity filter is slightly higher than two in both cases. This is consistent with experiments suggesting that two ions are permanently bound to the filter and arrival of a third ion triggers a fast conduction event.

View larger version (20K): [in this window] [in a new window]   FIGURE 9   Average number of potassium ions in a channel with an intrapore radius of 3.5 Å. The histograms are constructed in the same way as in Fig. 4. The distribution of ions in the channel are obtained in the presence of an applied field of +2 × 107 V/m (A) and 2 × 107 V/m (B).

Current-voltage relationships

In Fig. 10, we show how the current-voltage curves change with the intrapore radii, which are obtained using 300-mM solutions in both reservoirs. The curve obtained from the 4.5-Å channel is illustrated in Fig. 11 A. The curves derived from five different channel radii all reveal a distinct feature. The relationship is approximately linear when the applied potential is less than 100 mV, but it deviates systematically from Ohm's law with a further increase in the membrane potential. The degree of this nonlinearity increases systematically as the radius of the intracellular gate is reduced. This nonlinearity results from the presence of an energy barrier in the channel. Intuitively, a barrier is less of an impediment to an ion when the driving force is large. Rectification arises from the fact that ions moving into and out of the cell see different barriers. The solid lines fitted through the data points in Fig. 10 are simply meant to guide the eye. The slight shift of the I-V curves from the origin is due to the computational errors, which are quite substantial at low currents.

View larger version (11K): [in this window] [in a new window]   FIGURE 10   Current-voltage relationships. The current measured at various applied potentials is obtained with symmetrical solutions of 300 mM in both reservoirs. Each datum point in this and the following two figures is obtained from a simulation period of 3.2 µs. The intrapore radii of four model channels used are 3 ( in A), 3.5 ( in A), 4 ( in B), and 5 Å ( in B).

View larger version (11K): [in this window] [in a new window]   FIGURE 11   Current-voltage relationships obtained with the 4.5-Å channel with varying strengths of mouth dipoles. (A) Charge on each dipole guarding the intracellular entrance to obtain current-voltage relationship shown in open circles ± 0.4 e. The curve obtained with the charge on the dipole eliminated is shown in filled circles. In both curves, the dipoles near the channel entrance from the extracellular space carry the full unit charges. (B) The control curve (), reproduced from A, is compared with the curve obtained with the charge on the dipole guarding the extracellular entrance eliminated (). The charge on each intrapore mouth dipole is ±0.4 e, as in the control curve. Each data point for filled circles is derived from the simulation period of 3.2 µs.

In all the curves illustrated in Fig. 10, the outward current is approximately equal to the inward current at any given potential. This symmetry in the current-voltage relationship can be broken by changing the strengths of the dipoles guarding the channel entrances. To illustrate this point, we construct the current-voltage curves under various conditions using the 4.5-Å channel as an example. One such curve is shown in Fig. 11 A, where we remove the charges on the four dipoles located near the intracellular entrance. These charges represent the E₁₁₈ and R₁₁₇ residues in the detailed model and yield maximal conduction when set to q_i = 0.4 e. When the dipoles are eliminated, the outward current is severely attenuated (filled circles in Fig. 11 A). In other words, the channel under this condition is inwardly rectifying. For comparison, the control curve obtained in the presence of ±0.4 e on the dipoles is superimposed (open circles). In all previous simulations, each of the four dipoles placed near the entrance from the extracellular space, representing the D₈₀ and R₆₄ residues carries ±e. The current-voltage curve shown in Fig. 11 B is obtained with the charges on these dipoles removed, while keeping the charges of ±0.4 e on the intracellular dipoles (filled circles). For comparison, the control curve, shown in Fig. 11 A, is reproduced (open circles). Under this condition, the magnitude of outward currents remains unchanged, but the inward current is approximately halved, thus the channel becomes outwardly rectifying.

The shape of the current-voltage relationship revealed with these simulations is in disagreement with that obtained experimentally from the KcsA K⁺ channels inserted in planar lipid bilayer membranes (Schrempf et al., 1995; Heginbotham et al., 1999). The results of our simulations reveal that the current-voltage relationships are superlinear for the channels with small intrapore radii (Fig. 10 A) or when the charges on the mouth dipoles are reduced or eliminated (Fig. 11). The relationship is linear when the intrapore radius is large and the strengths of the mouth dipoles are optimum. In nearly all published data to date, with one exception (Meuser et al., 2001), the measured current is found to saturate at large applied potentials. However, it is suggested by Miller (1999) that such saturation invariably occurs due to "specific block by exogenous molecules," rather than due to "intrinsic ionic diffusion properties." In the absence of any blocking ions or molecules in the solution, we expect from theoretical considerations that the current-voltage relationship will deviate from Ohm's law, becoming superlinear. There are already some experimental indications for deviations from Ohm's law (for example, see Tyerman et al., 1992; Meuser et al., 2001). It will be of interest to test this prediction experimentally on a variety of other types of potassium channels, ensuring that ionic solutions contain no blocking agents and the applied potential is pushed beyond the usual range. If such deviations do occur, one may obtain an estimate of the barrier height in the channel using the formalism of Chung et al. (1999, Eq. 14). We note that the inward and outward rectification of the channel can be brought about by adjusting the ionization state of the charge residues near the channel gates. The current-voltage curves similar to those illustrated in Fig. 11 are observed in a number of ionic channels. In some cases (e.g., inwardly rectifying potassium channels), the rectification arises from block of the channel by internal Mg²⁺ or polyamines. Whether the rectification mechanism suggested here plays a role in other channels remains to be investigated.

Conductance-concentrations curves

Experimentally, the current I across the potassium channel first increases with an increasing ionic concentration [K] and then saturates, leading to a current-concentration relationship of the Michaelis-Menten form

(3)

so that the current approaches the saturation current I_max when [K] K_s, the half-saturation value. Theoretically, the conductance-concentration curve is expected to saturate if the transport of ions across the channel is determined by two independent processes, one of which depends on ionic concentrations in the reservoirs and one that does not. In the potassium channel, the time the energy well near the intracellular entrance of the channel remains vacant depends on ionic concentrations as well as the potential difference across the channel, whereas the time it takes for an ion to traverse the channel depends solely on the applied electric field (see Chung et al., 1999).

The magnitudes of current across the channel plotted against the concentrations of potassium ions in the reservoirs, as shown in Fig. 12, have the same shape as those observed experimentally. The five curves in Fig. 12 represent the outward current obtained from the channel with radii 3.0, 3.5, 4, 4.5, and 5 Å. The applied electric field used for all four curves is 2 × 10⁷ V/m. The conductance increases rapidly with an increasing ionic concentration initially and then saturates with a further increase in concentration. The half-saturation values of the curves are between 150 and 200 mM. The experimental K_s values for various types of potassium channels ranges from 40 mM for inward rectifiers (Stampe et al., 1998) to 300 mM for Shaker K⁺ channels (Heginbotham and MacKinnon, 1993).

View larger version (16K): [in this window] [in a new window]   FIGURE 12   Current-concentration curves. The outward current is measured with symmetric solutions of varying concentrations in the two reservoirs. The lines fitted through data points are calculated with Eq. 3. The curves obtained with the channels with five different radii are obtained with an applied field of 2 × 107 V/m.

	CONCLUSIONS

TOP ABSTRACT INTRODUCTION MATERIALS AND METHODS RESULTS AND DISCUSSION CONCLUSIONS REFERENCES

The family of potassium channels appear to have a common selectivity filter structure, but individually they exhibit a very diverse range of conductance properties. To understand this feature, we have constructed a simplified model of potassium channels starting from the crystal structure of KcsA and molding with MD to create an open state. In the simplified model, the protein is represented with a homogeneous, low-dielectric medium and the charge residues by various dipole groups. Both the electrostatic and conductance properties of the simplified channel are shown to be qualitatively similar to the atomic-detail model. We thus conclude that many of the salient properties of an ion channel can be deduced from a smoothed, simplified model, provided such a model captures the essential features of the real structure. Corollary to this statement is that, if the results of BD simulations closely reproduce experimental properties of a model channel, then the model used for simulations may have captured some of the structural features of the channel protein.

Having justified the use of a simplified model channel as a template, we then investigate the possible structural differences that could give rise to different potassium channels found in nature. We systematically change the radius of the intracellular pore entrance, leaving the dimensions of the selectivity filter and cavity unaltered. As the intrapore radius is increased from 2 to 5 Å, the channel conductance (at 245 mV) changes from 0.7 to 197 pS. By examining the energy profiles and the probabilities of ion occupancies in various segments of the channel, we deduce the rate-limiting step for conduction in the potassium channels. Ion distributions reveal that the selectivity filter is occupied by two K⁺ ions most of the time. A conduction event is triggered when a third ion climbs over the residual energy barrier located between the cavity and the intracellular mouth and moves toward the selectivity filter. This barrier is the rate-limiting step in the permeation process: as its height increases (with a decreasing intrapore radius), the channel conductance drops exponentially.

We have deduced some of the properties of the potassium channels that can be measured using patch-clamp techniques. These include the current-voltage curves and conductance-concentration relationships. The current-voltage curves are approximately linear when the applied potential is less than 100 mV, but for the channels with small intrapore radii, they deviate systematically from Ohm's law with a further increase in the membrane potential. When the charges on the mouth dipoles guarding the entrance of one side of the channel are reduced from the optimum values or eliminated altogether, asymmetries between the inward and outward currents result. The current increases with increasing ionic concentration and then saturates, leading to a current-concentration relationship of the Michaelis-Menten form. The half-saturation value is roughly the same for the channels with different radii but shifts upward when the driving force is increased.

We stress that the generic model used in this study is only a first step toward modeling of specific potassium channels. The intrapore shape and charges on it need to be further refined by exploiting the available mutation and electrophysiological data on target channels. For example, additional charged residues in the intrapore region are known to play an important role in determining the permeation properties of inward rectifiers (Thompson et al., 2000; Kubo and Murata, 2001). Reproducing detailed properties of specific potassium channels thus poses further challenges that will be taken up in future studies.