Mechanisms for the Ca2+-dependent gating of single large-conductance Ca2+-activated K+ channels from cultured rat skeletal muscle were developed using two-dimensional analysis of single-channel currents recorded with the patch clamp technique. To extract and display the essential kinetic information, the kinetic structure, from the single channel currents, adjacent open and closed intervals were binned as pairs and plotted as two-dimensional dwell-time distributions, and the excesses and deficits of the interval pairs over that expected for independent pairing were plotted as dependency plots. The basic features of the kinetic structure were generally the same among single large-conductance Ca2+-activated K+ channels, but channel-specific differences were readily apparent, suggesting heterogeneities in the gating. Simple gating schemes drawn from the Monod- Wyman-Changeux (MWC) model for allosteric proteins could approximate the basic features of the Ca2+ dependence of the kinetic structure. However, consistent differences between the observed and predicted dependency plots suggested that additional brief lifetime closed states not included in MWC-type models were involved in the gating. Adding these additional brief closed states to the MWC-type models, either beyond the activation pathway (secondary closed states) or within the activation pathway (intermediate closed states), improved the description of the Ca2+ dependence of the kinetic structure. Secondary closed states are consistent with the closing of secondary gates or channel block. Intermediate closed states are consistent with mechanisms in which the channel activates by passing through a series of intermediate conformations between the more stable open and closed states. It is the added secondary or intermediate closed states that give rise to the majority of the brief closings (flickers) in the gating.
Large-conductance Ca2+-activated K+ channels (maxi-K or BK channels),1 which are activated by micromolar concentrations of intracellular Ca2+ (Ca2+i) and by membrane potential, are present in a wide variety of tissues (reviewed by Rudy, 1988; Marty, 1989; McManus, 1991; Latorre, 1994). BK channels act to reduce membrane excitability by allowing K+ efflux though their opened pores to drive the membrane potential more negative. A large number of studies has led to the development of progressively more detailed models that can describe with increasing fidelity various features of the gating of BK channels (see McManus and Magleby, 1991; Wu et al., 1995; Cox et al., 1997, for recent models, and the references therein for preceding models). The minimal model of McManus and Magleby (1991), with three open and five closed states can account for many of the basic features of the Ca2+ dependence of the detailed single-channel gating kinetics over a 400-fold range of open probability (Popen = 0.00137–0.53). The 10-state minimal model of Wu et al. (1995) added two additional closed states to account for activity at high Ca2+i, and the 10- and 12-state minimal models of Cox et al. (1997) can account for the average Po (but not the single-channel kinetics) over a wide range of Ca2+i and voltage. The alpha subunits of BK channels assemble as a tetramer (Shen et al., 1994), and many of the above models that have been examined for the gating of BK channels can be related to allosteric mechanisms that have been proposed for tetrameric proteins, where each subunit can bind a single ligand. The allosteric model of Monod et al. (1965) has 10 states with concerted transitions, whereas the model of Koshland et al. (1966) has five states with sequential transitions. Both the concerted and sequential models are contained within the general 35-state model proposed by Eigen (1968) for tetrameric proteins.
The purpose of our present study is to investigate further the Ca2+ dependence of the detailed single-channel gating kinetics to examine which of the various allosteric models are consistent with the gating. As in our previous study (McManus and Magleby, 1991), simultaneous (global) maximum likelihood fitting of data obtained at several different Ca2+i was used to estimate the most likely rate constants and to rank the examined models, but with the fitting of two dimensional (2-D) rather than 1-D dwell-time distributions, to take the correlation information between adjacent open and closed interval durations into account (Fredkin et al., 1985; Magleby and Weiss, 1990b).
While full maximum likelihood fitting provides one of the best methods to rank models (Horn and Lange, 1983; Qin et al., 1997), it gives little detailed information about where the kinetic schemes may be inadequate, or insight into how they might be modified to better describe the data. To overcome these difficulties, the detailed kinetic information contained in the single-channel current record, the kinetic structure, was extracted and displayed as 2-D dwell-time distributions and dependency plots (Magleby and Song, 1992). Comparison of the observed kinetic structure to that predicted by the various kinetic schemes was then used to assess the ability of the top ranked models to account for the gating kinetics and to provide insight into how to modify the models to improve further the description of the data. Our findings identify minimal gating mechanisms for single BK channels that can account for both the Ca2+ dependence of the kinetics and the correlations between the durations of adjacent open and closed intervals. These models have additional brief closed states added to the model of McManus and Magleby (1991), which is a subset of the Monod-Wyman-Changeux model. These additional states can be located as intermediate closed states within, or as secondary closed states beyond, the activation pathway. Most of the flickers (brief closed intervals) are found to arise from these additional closed states. Our observations suggest that the Monod-Wyman-Changeux and Koshland-Nemethy-Filmer models are inadequate for BK channels, and that models based on either extensions of these models or on the more general model proposed by Eigen (1968) or extensions of Eigen's model may be more appropriate. A preliminary report of some of these findings has appeared (Rothberg and Magleby, 1998).
Materials And Methods
Currents flowing through single large-conductance Ca2+-activated K+ channels in patches of surface membrane excised from primary cultures of rat skeletal muscle (myotubes) were recorded using the patch clamp technique (Hamill et al., 1981). Pregnant rats were killed by CO2 inhalation, and cultures of myotubes were prepared from the fetal skeletal muscle as described previously (Barrett et al., 1982; Bello and Magleby, 1998). All recordings were made using the excised inside-out configuration in which the intracellular surface of the patch was exposed to the bathing solution. Kinetic analysis was restricted to patches containing a single BK channel. Single-channel patches were identified by observing openings to only a single open channel conductance level during several minutes of recording in which the open probability was >0.4. Experiments were performed at room temperature (22–24°C). All recordings were at +30 mV with the intracellular membrane surface positive.
The solutions bathing both sides of the membrane contained either 150 mM KCl and 5 mM TES (N-tris(hydroxymethyl)methyl-2-aminoethane sulfonate) pH buffer, with the pH of the solutions adjusted to 7.0 (channels B06, B12, B14), or 144 mM KCl, 2 mM TES, and 1 mM EGTA, with the pH of the solutions adjusted to either 7.2 (channels M09 and M25) or 7.0 (M24). The solution at the intracellular side of the membrane also contained added Ca2+ (as CaCl2), to bring the free Ca2+ concentration at the intracellular surface (Ca2+i) to the indicated levels. No Ca2+ was added to the extracellular (pipette) solution. Ca2+i in the absence of EGTA was determined by atomic absorption spectrometry, and in the presence of EGTA was determined with a Ca2+-sensitive electrode as detailed previously (McManus and Magleby, 1988). Solutions were changed though the use of a microchamber (Barrett et al., 1982).
Recording and Measuring Interval Durations and Identifying Normal Activity
Single-channel currents were recorded on either FM tape (DC-20 kHz) or on a digital data recorder (DC-37 kHz). Data were then low-pass filtered with a four-pole Bessel filter to give a final effective filtering of typically 6–10 kHz (−3 dB) and sampled by computer at a rate of 80–200 kHz. The filtering gave dead times (the duration of an interval that would just reach the 50% current level) for the various channels of (μs): 28.5 B06, 22.9 B12, 17.9 B14, 36.0 M25, 28.6 M09, 24.6 M24. The sampled currents were then analyzed using custom programs written in the laboratory. The methods used to set the level of filtering to exclude false events that could arise from noise, measure interval durations with half-amplitude threshold analysis, test for stability, and identify modes using stability plots, have been described previously, including the precautions taken to prevent artifacts in the analysis (McManus et al., 1987; McManus and Magleby, 1988, 1989; Magleby, 1992). We thank Owen McManus (Merck Sharp and Dohme Research Laboratories, Manukau City, New Zealand) for data from channels M09, M24, and M25, which were originally obtained for the study by McManus and Magleby (1991). The analysis in the present study was restricted to channel activity in the normal mode, which typically involves 96% of the detected intervals (McManus and Magleby, 1988). Activity in modes other than normal, including the low activity mode (Rothberg et al., 1996), was removed before analysis. Since the analysis in this study takes into account the correlation information in the durations of adjacent open and closed intervals, the sites of any removal of intervals due to activity in modes other than normal or artifacts associated with transitions to subconductance levels, were marked to avoid the later juxtaposition of open and closed intervals that were not adjacent in the original record.
Log Binning and Plotting of 2-D Dwell-Time Distributions
Two types of 2-D dwell-time distributions were generated. The first was 2-D frequency histograms of each pair of successive (adjacent) open and closed intervals. These distributions were used for the maximum likelihood fitting to determine the minimal numbers of components (states) and also for evaluating kinetic schemes to estimate the most likely rate constants and obtain likelihood estimates. The second type was surface plots constructed using interpolated smoothing of the histograms. The 2-D distributions presented in the figures are surface plots.
The first step in generating a 2-D frequency histogram for the dwell-time distribution was to bin adjacent open and closed intervals. Every open interval and its following (adjacent) closed interval were binned and every closed interval and its following (adjacent) open interval were also binned, with the logs of the open and closed interval durations of each pair locating the bin on the x and y axes, respectively. Each interval was thus binned twice, but with a different adjacent interval. Including open–closed and closed–open interval pairs in each distribution assumes microscopic reversibility, an assumption that appears consistent with the data (Song and Magleby, 1994). The 2-D frequency histograms for the 2-D dwell-time distributions were binned at a resolution of 10 per log unit. Further details on log binning of 2-D dwell-time distributions may be found in Magleby and Weiss (1990b) and Rothberg et al. (1997).
The surface plots for display of the 2-D dwell-time distributions were constructed from the 2-D frequency histograms in a series of steps. The first step was to smooth the histograms using a 2-D moving bin average with three bins per side, with the number of events in each bin weighted as the inverse of the distance from the central bin. Thus, the numbers of events in the four corner bins in the three-by-three moving array were multiplied by 0.707 before being added to the events in the other bins of the moving array. The total was then divided by 7.828 (4 × 0.707 to weight the corner bins plus 5 × 1 to weight the center and noncorner bins) to obtain the weighted average for the bin in the position of the center bin in the new smoothed distribution. The process was then repeated for all bins in the unsmoothed distribution to obtain the values for the new smoothed distribution. The Sigworth and Sine (1987) transform, which plots the square-root of the numbers of observations per log bin, where the bin widths are constant on a log scale, was then applied to the smoothed distribution.
Once the data were transformed, the 2-D surface plots for display were generated with the program Surfer (Golden Software, Golden, Colorado). The interpolation for the gridding with Surfer was performed using the inverse distance to a power method with smoothing, where the power was 2.0 and the smoothing factor was 0.1. Applying these smoothing procedures to distributions generated from different numbers of simulated intervals indicated that the smoothing procedures reduced features that might be expected to arise from stochastic variation, while having little effect on the basic features. The smoothing procedures were used only for visual display. The fitting was performed on the 2-D frequency histograms without averaging or smoothing. To simplify the writing, the text will refer to the fitting of 2-D dwell-time distributions presented in the figures, when, in reality, it was the 2-D frequency histograms that were fitted.
With filtering, detected intervals with durations less than approximately twice the dead time are narrowed (McManus et al., 1987; Colquhoun and Sigworth, 1995). For the fitting of kinetic models using 2-D frequency histograms, the measured durations of these intervals were corrected to their estimated true durations before binning and fitting, using the numerical method in Colquhoun and Sigworth (1995). For the surface plots presented in the figures, the measured durations were not corrected for narrowing before binning and plotting.
Dependency plots were constructed from the 2-D dwell-time distributions as detailed in Magleby and Song (1992). Briefly, the dependency for each bin of open-closed interval pairs with mean durations tO and tC is
where NObs(tO,tC) is the experimentally observed number of interval pairs in bin (tO,tC), and NInd(tO,tC) is the calculated number of interval pairs in bin (tO,tC) if adjacent open and closed intervals pair independently (at random). The method of calculating expected frequencies for criteria that are independent (contingency tables) is a common statistical procedure (see Mendenhall et al., 1981). The expected number of interval pairs in bin (tO,tC) for independent pairing is
where P(tO) is the probability of an open interval falling in the row of bins with a mean open duration of tO, and P(tC) is the probability of a closed interval falling in the column of bins with a mean closed duration of tC. P(tO) is given by the number of open intervals in row tO divided by the total number of open intervals in all rows, and P(tC) is given by the number of closed intervals in the column in tC divided by the total number of closed intervals in all columns. Since open and closed intervals are paired, the total number of open intervals is equal to the total number of closed intervals, which is equal to the total number of interval pairs in the 2-D dwell-time distribution.
Estimating the Most Likely Rate Constants for Kinetic Schemes
The most likely rate constants for the examined kinetic schemes were determined from the simultaneous fitting of the 2-D frequency histograms (dwell-time distributions) obtained at three different Ca2+i using an iterative maximum likelihood fitting procedure similar to the one detailed in McManus and Magleby (1991), except that 2-D dwell-time distributions replaced the 1-D dwell-time distributions, and the correction method of Crouzy and Sigworth (1990) for missed events due to filtering replaced our previous correction method. The steps in the fitting were: (a) for the given kinetic scheme and starting rate constants, an equivalent uncoupled kinetic scheme (Kienker, 1989) with additional kinetic states to account for missed events was calculated based on the dead time and Ca2+i; (b) the time constants and volumes of the 2-D components underlying the predicted 2-D dwell-time distributions for the given kinetic scheme and filtering were calculated from the equivalent kinetic scheme using 2-D Q-matrix methods (Fredkin et al., 1985; Colquhoun and Hawkes, 1995b); (c) the likelihood that the interval pairs in the observed 2-D dwell-time distribution were drawn from the predicted distribution was then calculated using the predicted underlying 2-D components, as detailed in Rothberg et al. (1997); (d) steps a–c were repeated for the 2-D distribution obtained at each different Ca2+i, and the global log likelihood for the simultaneously fitted 2-D dwell-time distributions was then the sum of the log likelihoods for the individual distributions; and (e) the rate constants were then changed using a maximization routine. Steps a–e were repeated until the rate constants for the given scheme and dead time were found that maximized the likelihood.
Precautions were taken during the fitting to diminish the chance that the rate constants for a given fitted scheme represented a local maximum on the likelihood surface. For example, schemes were typically refit using different initial rate constants, and the size of the step change for each rate constant was varied and periodically reset during the maximization routine to increase the possibility of jumping over local maxima. In spite of these precautions, we cannot exclude that more likely fits might be found in some cases.
Estimating the Significance of the Dependencies
The significance of the dependencies was obtained by comparing the numbers of interval pairs in the various bins of the observed 2-D dwell-time distribution with the number expected if adjacent open and closed intervals paired independently. The comparison was made using a moving paired t test for nine bins at a time in corresponding three-by-three arrays from the observed and expected distributions. After each comparison, both arrays were moved one bin, until the entire surface of the 2-D distribution was covered. The calculated P value was determined from a t table with eight degrees of freedom, and then converted to the log of the P value times the sign of the dependency. This dependency significance was then plotted at the centers of the moving three-by-three arrays to generate 2-D dependency significance plots. With this transform, dependency significance values >1.3 or <−1.3 would indicate P values <0.05. Heavy lines at ±1.3 were included on the dependency significance plots to indicate when the dependencies were significant for P < 0.05.
Estimating the Theoretical Best Description of the 2-D Dwell-Time Distributions
To evaluate models, it is useful to have an estimate of the theoretical best description of the dwell-time distributions. This theoretical best description can then be compared with the best description generated by any given kinetic scheme in order to evaluate how well the kinetic scheme describes the data. If the assumption is made that the gating of the BK channel is consistent with a discrete state Markov process, such that the rate constants do not change with time (McManus and Magleby, 1989; Petracchi et al., 1991), then two different methods can be used to obtain an estimate of the theoretical best description of the 2-D dwell-time distributions that would be obtained if the discrete state Makov gating mechanism were known.
In the first method, the 2-D dwell-time distributions were fitted with sums of 2-D exponential components with all free parameters, except for the volume of one component, as the volumes of the components must sum to 1.0 (Rothberg et al., 1997). The number of components was increased until there was no longer a significant increase in likelihood. The maximum likelihood for this fit would then approximate that of the theoretical best description for a discrete state Markov model fit to the exact same data.
In the second method, an uncoupled kinetic scheme equivalent to the unknown gating mechanism was used to estimate the theoretical best fit to the data. This approach is based on an extension of the observation of Kienker (1989), who found that any given kinetic scheme can be transformed into an equivalent uncoupled kinetic scheme. Since the form of the uncoupled scheme depends only on the number of open and closed states, then it follows that the uncoupled scheme for a channel can be determined without knowing the gating mechanism, provided that the numbers of open and closed states are known. Although the gating mechanism of the uncoupled scheme is different from the unknown gating mechanism, the uncoupled scheme with appropriate rate constants should give descriptions of the single-channel data that are identical to those that would be obtained from the (unknown) underlying kinetic scheme. Hence, fitting a 2-D dwell-time distribution with an uncoupled scheme should give the same theoretical best description of the distribution as the unknown kinetic scheme, assuming a discrete state Markov model and provided that both schemes have the same number of states. To estimate the theoretical best likelihood for the simultaneous fitting of 2-D dwell-time distributions obtained under different experimental conditions, each distribution was fitted separately with an uncoupled scheme, and then the log likelihoods for the separate distributions were summed together.
Estimating the theoretical best likelihood by fitting with uncoupled schemes has an advantage over fitting with sums of 2-D components in that uncoupled schemes can be used to simulate single-channel data with filtering and noise, provided that none of the rate constants in the fitted uncoupled schemes are negative, which appears to be the case so far. While the uncoupled schemes can give an estimate of the theoretical best likelihood, they do not have predictive value beyond the specific experimental conditions for the data they are fitted to, as there are no Ca2+- or voltage-dependent rate constants in the uncoupled schemes.
Ranking the Kinetic Schemes
Normalized likelihood ratios (NLR) have been used to indicate how well any given kinetic scheme describes the 2-D dwell-time distributions when compared with the theoretical best description of the data. Normalization accounts for the differences in numbers of interval pairs among experiments, so that comparisons can be made between channels. The normalized likelihood ratio per 1,000 interval pairs is defined as
where ln S is the natural logarithm of the maximum likelihood estimate for the observed 2-D dwell-time distributions given the kinetic scheme, ln T is the natural logarithm of the maximum-likelihood estimate for the theoretical best description of the observed distributions, and N is the total number of fitted interval pairs (events) in the observed dwell-time distributions (McManus and Magleby, 1991; Weiss and Magleby, 1992).
A value of 1.0 for the NLR1000 indicates that the given kinetic scheme describes the observed distributions as well as theoretically possible for a discrete state Markov model. A value of 0.05 would indicate that the probability that the observed data were generated by the given kinetic scheme is only 5% per 1,000 interval pairs of the probability that the observed data were derived from the theoretical best description of the distributions.
The NLR gives a measure of how well different kinetic schemes describe the distributions, but it cannot be used to directly rank schemes, since no penalty is applied for the numbers of free parameters. To overcome this difficulty, the Schwartz criterion has been used to apply penalties and rank models (Schwarz, 1978; Ball and Sansom, 1989). The Schwarz criterion (SC) was calculated from
where L is the log-likelihood value, F is the number of free parameters, and N is the number of interval pairs. The scheme with the smallest SC is the top ranked scheme.
Experimental single-channel data is distorted by the combined effects of noise and low-pass filtering. Thus, to make valid comparisons between the observed distributions and the distributions predicted by the kinetic models, simulated single-channel current records were generated with filtering equivalent to that used to analyze the single-channel current and with noise similar to that in the single-channel current. The simulated single-channel currents were then analyzed in the same manner used to analyze the experimental currents. The method used to simulate single-channel currents with true filtering and noise is detailed in Magleby and Weiss (1990a).
Currents flowing through a single large-conductance Ca2+-activated channel in an inside-out patch of membrane excised from a cultured rat skeletal muscle cell are shown in Fig. 1, A and B, at two different time bases. The complexity of the underlying gating process is reflected in the wide range of the durations of the open and closed intervals and the apparent grouping of the intervals into bursts. Because of the stochastic nature of single-channel gating (Colquhoun and Hawkes, 1995a), extracting the essential kinetic information about the underlying gating process requires large amounts of stable single-channel data to overcome the stochastic variation. Fig. 1,C and D, presents stability plots of the mean open and closed interval durations during activity in the normal mode, which includes ∼96% of the intervals (McManus and Magleby, 1988). The stability plots shown in Fig. 1 are based on 57.6 s of stable data after artifacts and transitions to modes other than normal were removed. These stability plots indicate that the analyzed data are reasonably stable, and are representative of the data analyzed in this study to investigate Ca2+-dependent gating mechanisms.
2-D Dwell-Time Distributions
For channels that gate between two conductance levels, open and closed, 2-D dwell-time distributions contain essential kinetic information from the single-channel current record, including correlation information that gives information about transition pathways among states (Fredkin et al., 1985; Rothberg et al., 1997). Fig. 2 shows 2-D dwell-time distributions for six single BK channels, each from a different inside-out patch of surface membrane. The membrane potential in each case was +30 mV and the Ca2+i was selected to give a Po near 0.5. The 2-D dwell-time distributions plot how frequently open intervals of a specified duration occur next (adjacent) to closed intervals of a specified duration. The log of the durations of each adjacent open and closed interval locate the bin on the x and y axes, respectively, and the z axis plots the square root of the number of observations per bin (see materials and methods).
The 2-D dwell-time distributions in Fig. 2 can be described by the sums of 2-D exponential components, where the number of 2-D components is given by the product of the numbers of open and closed states (Fredkin et al., 1985; Rothberg et al., 1997). Since BK channels typically enter a minimum of three to four open and five to seven closed states during normal activity (McManus and Magleby, 1988), there would be from 15–28 possible 2-D components underlying each 2-D dwell-time distribution. The square-root transformation (Sigworth and Sine, 1987) used for the 2-D dwell-time distributions would generate peaks at the time constants of the 2-D exponential components (Rothberg et al., 1997). However, only the components with the largest volumes or whose time constants are well separated from the other components would generate visually detectable peaks.
To facilitate reference to the various peaks and regions of the 2-D plots, the 2-D distributions in Fig. 2 and in the subsequent figures are divided into six general regions indicated by the numbers 1–6 in the figures and referred to as #1–#6 in the text. For example, #1, #2, and #3 indicate the regions of brief openings adjacent to brief closings, intermediate closings, and long closings, respectively, while #4, #5, and #6 indicate the regions of long openings adjacent to brief closings, intermediate closings, and long closings, respectively.
The highest peak (#4) in the 2-D dwell-time distributions in Fig. 2 is located in the same general position for each channel and indicates that the most frequent interval pairs for all the examined channels consisted of long (∼2-ms) openings adjacent to brief (∼0.05-ms) closings. These dominant interval pairs are readily apparent in the single-channel record in Fig. 1 B as longer open intervals adjacent to the brief closed intervals (flickers). Other peaks and inflections are also apparent in the plots, indicating that additional components can be detected visually. For example, each plot contains a visible peak indicating a component of long (∼2-ms) openings adjacent to long (∼10-ms) closings (#6), and a component of brief openings (∼0.1 ms) also adjacent to long closings (10 ms) (#3).
The 2-D dwell-time distributions in Fig. 2 indicate the relative frequency of occurrence of the various classes of adjacent open- and closed-interval durations that must be accounted for by kinetic gating mechanisms. The channels for Fig. 2 were selected to be representative of the more than 12 channels examined in this manner. Channels M25, M09, and M24 are channels 1, 2, and 5, respectively, in McManus and Magleby (1991), and were included to allow comparison of the 2-D analysis in this present study with the 1-D methods used previously.
Kinetic Similarities and Heterogeneities for BK Channels from the Same Preparation
While there are a number of basic similarities in the 2-D dwell-time distributions from the six different BK channels in Fig. 2, there are also a number of differences. For example, channels M25 and M09 have a prominent middle ridge (#5), indicating a component of long openings (∼2 ms) adjacent to intermediate duration closed intervals (∼0.5 ms). This component is less apparent or appears to be missing for the other four channels. Although there were some differences in the level of filtering and Po among the different channels (see Fig. 2, legend), it is unlikely that this would account for the differences in the 2-D dwell-time distributions as there was no evident relationship between the observed differences and the small differences Po, or the levels of filtering for the various channels.
The obvious kinetic differences among the channels in Fig. 2 are consistent with previous studies showing differences in Ca2+ sensitivity and/or gating among different native BK channels from the same tissue (McManus and Magleby, 1991; Wu et al., 1996). Since the six different channels in Fig. 2 were obtained from native tissue, it is possible that the differences in kinetics might reflect channels with different splice variations (Atkinson et al., 1991; Adelman et al., 1992; Butler et al., 1993; Lagrutta et al., 1994). Alternatively, other factors may be involved (see discussion) since expressed cloned channels without the potential for alternative splicing can also display kinetic differences among channels (Silberberg et al., 1996). Kinetic heterogeneity has been observed for other types of channels as well (e.g., Auerbach and Lingle, 1986; Patlak et al., 1986).
Displaying the Correlations between Adjacent Open and Closed Intervals with Dependency Plots
Although information on the correlation between adjacent open and closed intervals is contained within the 2-D dwell-time distributions, this information is not readily apparent from visual inspection. Dependency plots provide a means to extract this correlation information in a form that can give insight into the connections among open and closed states involved in the gating (Magleby and Song, 1992).
Dependency plots for the six channels shown in Fig. 2 are presented in Fig. 3. The plots present the fractional differences between the observed number of adjacent open and closed intervals of indicated durations and the hypothetical number that would be observed if all the open and closed intervals paired independently. Dependencies of +0.5 or −0.5 would indicate a 50% excess or deficit, respectively, of interval pairs over that expected for random pairing (see Eq. 1). Positive dependencies suggest that the open and closed states underlying the interval pairs in excess are effectively connected, and negative dependencies suggest that the open and closed states underlying the interval pairs in deficit are not effectively connected (Magleby and Song, 1992).
The dependency plots in Fig. 3 show some common kinetic features for the six different BK channels: a deficit of brief open intervals adjacent to brief closed intervals (#1), an excess of brief open intervals adjacent to long closed intervals (#3), an excess of long open intervals adjacent to brief closed intervals (#4), and a deficit of long open intervals adjacent to long closed intervals (#6). It will be shown in the next section that these specific dependencies are significant. Thus, the basic features of the dependency plots in Fig. 3 suggest that kinetic models for the gating should contain dominant transition pathways between the open and closed states underlying the brief open intervals and the long closed intervals (#3), and between the open and closed states underlying the long open intervals and brief closed intervals (#4). These dominant transition pathways would generate the positive dependencies. In addition, there should not be dominant transition pathways between the open and closed states underlying the brief open and closed intervals (#1) and between the open and closed states underlying the long open intervals and long closed intervals (#6). Whether a transition pathway is dominant or not depends on the relative probability of whether that pathway is taken among the possible pathways from any given state.
Interestingly, the dependency plot of channel M09 showed features that were not observed in any of the other BK channels: most notably, a smaller excess of brief open intervals adjacent to longer closed intervals (#3). This atypical kinetic structure of channel M09 is consistent with differences in the gating mechanism for this channel when compared with four other channels, determined in a previous study (channel 2 vs. channels 1, 3, 4, and 5 in McManus and Magleby, 1991). That channel M09 is atypical is readily apparent from the dependency plots in Fig. 3 obtained at a single Ca2+i, indicating the power of dependency plots. The previous determination that channel M09 was atypical required hundreds of hours of analysis of 1-D dwell-time distributions obtained at three or more Ca2+i for each channel.
Determining the Significant Features of Dependency Plots
While the basic features of the dependency plots were consistent among most channels, there were also channel-specific features in the plots. Therefore, it was of interest to determine which features were part of the kinetic structure and which might have arisen from factors such as stochastic variation, noise in the single-channel records, and distortions produced by low-pass filtering. The significance of the dependencies were estimated by two different approaches: using simulation and applying a paired t test.
The first approach used simulation to estimate the magnitude of the expected variations in the dependency that would arise from stochastic variation, filtering, and noise. A 2-D dwell-time distribution and associated dependency plot were simulated for Scheme I, a gating mechanism that would give theoretical dependencies of zero (Magleby and Song, 1992). Scheme I was first fitted to the 2-D dwell-time distribution for channel B06 to estimate the most likely rate constants. These rate constants were then used with Scheme I to simulate a single-channel current record with noise and filtering like that in the experimental data and with the same number of intervals as for channel B06. The simulated current record was then analyzed to obtain the 2-D dwell-time distribution and dependency plot shown in Fig. 4. The deviations of the dependency plot from zero in Fig. 4 B then give an estimate of the variations that would be expected due to the combined effects of noise, filtering, and stochastic variation, since the expected theoretical dependencies for Scheme I would be zero.
From Fig. 4 and similar simulations of this type, the magnitudes of the variations of the dependency from zero were found to depend on the location in the plot. There was little deviation from zero for the dependency of long open intervals adjacent to brief closed intervals (#4) because of the large numbers of interval pairs that contributed to the calculation of dependency for this location, as seen in the 2-D dwell-time distribution. Elsewhere in the plot, the dependencies that would be expected to arise from stochastic variation typically fell within ±0.2 from zero.
A comparison of the predicted dependency plot in Fig. 4,B to the observed dependency plots in Fig. 3 indicates that Scheme I is inconsistent with the gating of BK channels. Nevertheless, the 2-D dwell-time distribution predicted by Scheme I in Fig. 4,A appears similar (but not identical) to the observed dwell-time distribution in Fig. 3 for channel B06. Hence, the ability of a model to approximate the 2-D dwell-time distribution of the data does not necessarily establish that even the basic features of the proposed gating mechanism are correct. 1-D distributions can be even less sensitive for model discrimination (Magleby and Weiss, 1990b).
The second approach to estimate the significance of the dependencies involved a direct calculation of significance. Figs. 4,C and 5 plot the statistical significance of the dependencies. The significance was estimated by comparing the numbers of intervals in the observed 2-D dwell-time distribution with the number expected if adjacent open and closed intervals paired independently of one another. The comparison was made using a paired t test (details in materials and methods). The distributions of dependency significance in Fig. 5 plot the significance of the dependencies in Fig. 3 as the logarithm of the estimated P value, which is then multiplied by the sign of the dependency to indicate whether the paired intervals are in excess or deficit. The heavy lines on the plots at −1.3 and 1.3 indicate a significance level of P = 0.05. Absolute values of dependency significance >1.3, 2, 3, and 4 would indicate P < 0.05, < 0.01, < 0.001, and < 0.0001, respectively.
The dependency significance plots in Fig. 5,A are for the same six channels and orientation as in Fig. 3 (front view). Fig. 5,B presents the backside views of the dependency significance plots for two of the six channels to show the significant deficit of long openings adjacent to long closings (#6). Similar significant deficits were seen for the other four channels. It is important not to confuse the significance of the dependency with the magnitude of the dependency. Fig. 3 shows the magnitude of the dependencies. Fig. 5 shows whether the indicated magnitudes are significant. Fig. 4 C provides an independent measure of the applied significance test, showing that none of the dependencies arising from stochastic variation were significant, as would be expected for Scheme I.
Kinetic Structure of BK Channels
The 2-D dwell-time distributions and dependency plots in Figs. 2 and 3 and the significance of the dependencies in Fig. 5 are representative of dependency plots obtained from more than 12 channels. These plots present the essential kinetic information contained in the single-channel current records, indicating the kinetic structure of the BK channels. It is this information that must be accounted for by proposed gating mechanisms.
Idealized Dependency Plots from Single-Channel Data
It would be useful if there were a means to eliminate the variation in dependency plots arising from the analysis of limited amounts of single-channel data. We have developed an approximate means to do this, by fitting the data with an uncoupled (generic) kinetic scheme. Since the uncoupled scheme allows direct transitions from each open state to each closed state (Kienker, 1989), the correlations between the adjacent open and closed intervals can be described by such a scheme. The uncoupled kinetic scheme for any discrete state Markov model with four open and six closed states is given by Scheme II. The scheme is uncoupled because there are no direct transition pathways from one open state to another or from one closed state to another.
Although the actual kinetic scheme for channel B06 is not known, if the data are described by four open and six closed exponential components, which is the case for the distribution in Fig. 3, then Scheme II with most likely rate constants should give the same description of the experimental data as the unknown kinetic scheme (see materials and methods), and should thus be capable of describing the kinetic structure. Scheme II can then be used to generate simulated single-channel currents with different levels of noise and filtering, to examine the effects of these variables on the kinetic structure.
Scheme II was fitted to the 2-D dwell-time distribution for channel B06 to obtain the most likely rate constants. Scheme II was then used to simulate single-channel currents with filtering and noise similar to that in the experimental data, and also without filtering and noise. 1,000,000 detected intervals were simulated in each case to reduce stochastic variation to negligible levels. The simulated single-channel currents were then analyzed to generate the idealized 2-D dwell-time distributions and dependency plots presented in Fig. 6. These idealized distributions give an estimate of what the experimental distributions would look like without stochastic variation (Fig. 6,A), and without filtering, noise, or stochastic variation (B). A comparison of the idealized distributions in Fig. 6,A to those for channel B06 in Figs. 2 and 3 suggest that the minor variations in the experimental dependency plots most likely arise from stochastic variation due to the analysis of limited amounts of data. A comparison of the distributions in Fig. 6,A with filtering and noise to those in Fig. 6 B without filtering and noise indicates that filtering and noise do not change the basic features of the kinetic structure, except for those features involving adjacent intervals in which one or both intervals have durations less than two dead times, where the filtering attenuates the durations (McManus et al., 1987; Colquhoun and Sigworth, 1995).
Calcium Dependence of the Kinetic Structure
McManus and Magleby (1991) have previously detailed the Ca2+ dependence of the 1-D dwell-time distributions of open and closed interval durations for BK channels from rat skeletal muscle. The Ca2+ dependence of the 1-D distributions for the seven additional BK channels analyzed in this manner for the present study (data not shown) were consistent with those from the BK channels analyzed previously. Increasing [Ca2+]i increased Po by increasing the mean open interval duration and decreasing the mean closed interval duration. The Ca2+-dependent shifts in mean interval durations arose mainly from decreases in the time constants and areas of the longer closed components and increases in the time constants and areas of the longer open components. In contrast to the changes in the time constants of the longer open and closed components, the time constants of the shortest open and closed components appeared relatively independent of Ca2+i.
To gain further insight into the Ca2+ dependence of the gating, the Ca2+ dependence of the kinetic structure was examined. Results are presented in Fig. 7 for a representative channel (B06) at three different Ca2+i of 5.5, 8.3, and 12.3 μM, which resulted in open probabilities of 0.061, 0.202, and 0.504. Examples of single-channel current records at each level of activity are shown in Fig. 7,A and the kinetic structures are shown in Fig. 7,B. The greater variation in the dependency at 5.5 and 8.3 μM Ca2+i is most likely due to the fewer intervals obtained for analysis at these distributions due to the lower levels of activity. The idealized kinetic structures obtained after removing the effects of stochastic variation (as discussed for Fig. 6,A) are shown in Fig. 7 C. These idealized plots show the dominant features of the kinetic structure.
The 2-D dwell-time distributions show a shift in the time constants of the longer closed intervals towards shorter durations with increasing Ca2+i, and also a shift in the frequency of occurrence of the longer closed intervals towards briefer closed intervals, and of briefer open intervals towards longer open intervals. The characteristic saddle-like appearance of the dependency plots, which indicates an inverse relationship between the durations of adjacent open and closed intervals (McManus et al., 1985), was maintained at each level of activity, suggesting that the basic underlying gating mechanism remained unchanged over the examined range of Po.
A Simple Gating Mechanism Approximates the Basic Features of the Ca2+ Dependence of the Kinetic Structure
The basic features of the Ca2+ dependence of the 1-D dwell-time distributions from BK channels from cultured rat skeletal muscle can be described by Scheme III, which contains three open and five closed states (McManus and Magleby, 1991). To test whether this scheme might also account for the basic features of the kinetic structure, 2-D distributions obtained at three different Ca2+i from each channel were simultaneously fitted to Scheme III to estimate the most likely rate constants for each channel. These rate constants were then used with Scheme III to simulate a single-channel current record that was then analyzed to obtain the predicted kinetic structure. The current record was simulated with filtering and noise like that in the experimental data, and 106 simulated intervals were analyzed for each distribution.
Fig. 8 shows the 2-D distributions and dependency plots predicted by Scheme III for channel B06. Comparison of these predicted kinetic structures to the observed and idealized kinetic structures in Fig. 7, B and C, shows that Scheme III approximates the basic features of the Ca2+ dependence of the kinetics. However, there are some clear differences between the observed and predicted distributions. Scheme III predicted a greater deficit of brief open intervals adjacent to brief closed intervals for 5.5 μM Ca2+i than was observed in the experimental data (#1). That is, Scheme III generated an insufficient number of brief open intervals adjacent to brief closed intervals. This underprediction by Scheme III prompted a search for adjacent brief open and closed intervals in the single-channel current record. Fig. 9 shows examples of such pairings during normal gating with 5.5 μM Ca2+i. Approximately 20– 30% of adjacent brief open and closed intervals were found at the beginnings and endings of bursts and 70– 80% were found within bursts.
In addition, Scheme III predicted a greater excess of brief open intervals adjacent to intermediate duration closed intervals than were observed, which was apparent at all three levels of Ca2+i in both the dependency plots (#2) and the 2-D dwell-time distributions. Finally, Scheme SIII underpredicted the observed excess of long open intervals adjacent to the brief closed intervals (#4).
Similar results were found for four additional channels examined in detail. Scheme III captured the basic features of the Ca2+ dependence of the kinetic structure while giving the same types of over- and underpredictions, often with greater differences than those detailed for channel B06 above.
The Dependency Plots Suggest how Scheme III Might Be Modified to Better Describe the Kinetic Structure
To gain possible insight into why Scheme III did not account for all the features of the kinetic structure, the mean lifetimes of the kinetic states in Scheme SIII were calculated from the most likely rate constants. The results are presented in Scheme III(8.3) for channel B06 with 8.3 μM Ca2+. The solid line encloses the major gating pathways, and the mean lifetimes of the various states are indicated (in milliseconds). For kinetic schemes with compound open and closed states, it can be difficult, if not impossible, to designate which states contribute to the observed components of interval durations (Colquhoun and Hawkes, 1981, 1995). Nevertheless, for certain schemes and rate constants, such as Scheme III, such assignments can be tentatively made for the purposes of investigating why the scheme did not account for the complete features of the kinetic structure. The assignments were made by changing the lifetimes of the states one at a time, in small amounts, to determine which components were affected in the calculated distributions. The details of this approach will be presented elsewhere.
Scheme III predicted too few brief open intervals adjacent to brief closed intervals in low Ca2+i (#1 in Figs. 7 and 8). For Scheme III, brief openings, which are mainly from O3, would occur infrequently adjacent to brief closings, which are mainly from C5, since transitions from O3 and C5 must pass through either the intermediate closed state, C6, or the intermediate open state, O2. Transitions through either of these intermediate states would extend the mean duration of the open or closed intervals so that the intervals would no longer be brief. To compensate for the inability of Scheme III to generate a sufficient number of brief open intervals adjacent to brief closed intervals, a transition pathway to an additional brief duration closed state could be added to O3, as in Scheme IIIA. This would increase the number of brief openings adjacent to brief closings.
Scheme III also generated an excess of brief openings adjacent to intermediate-duration closings (#2 in Figs. 7 and 8). This excess is likely to arise from direct transition between O3 and C6. While the lifetime of C6 is ∼0.4 ms, transitions such as O3-C6-C5-C6-O3 would double the mean duration of the closed interval to ∼1 ms. Transitions such as O3-C6-O3-C6 could also increase the observed mean duration of the closed interval associated with state C6, since some of the dwell-times in O3 would be too brief to detect due to the filtering, making the apparent closed interval longer due to missed events (Blatz and Magleby, 1986; Colquhoun and Sigworth, 1995). The brief duration closed state C11 connected to O3 in Scheme IIIA would increase the number of brief openings adjacent to brief closings. C11 would also act to decrease the excess of brief openings adjacent to intermediate closings (mainly from C6) by diverting some of the transitions from O3 away from C6.
Finally, Scheme III generated an insufficient number of brief closed intervals adjacent to long open intervals (#4 in Figs. 7 and 8). In Scheme III, brief closings adjacent to long openings would arise mainly from transitions such as O1-O2-C5-O2-O1. Adding the brief duration closed states C9 and C10 to O1 and O2, as in Scheme IIIA, would increase the number of brief closings adjacent to long openings, as was observed in the experimental data.
Thus, a comparison of the differences between the observed dependency plots and those predicted by Scheme III suggest that there may be additional brief closed states directly connected to the open states, as proposed in Scheme IIIA. Transitions to these brief closed states would add flickers to the single-channel current record. Since these additional closed states are not on the activation pathway, they will be referred to as secondary closed states, giving rise to secondary flickers. The primary flickers would involve transitions to C5 on the activation pathway. The addition of the secondary closed states would represent a minimal extension of Scheme III, so that the modified Scheme IIIA should still capture the basic features of the gating, while reducing some of the differences between the observed and predicted dependency plots. Further reason for investigating closed states beyond the activation pathway comes from the observation of Wu et al. (1995) that such states were required to account for activity with high Ca2+i for the gating for BK channels in hair cells. Such secondary closed states could arise from either a secondary gate or channel block, as will be considered in the discussion.
Scheme IIIA, which has Secondary Closed States, Improves the Description of the Kinetic Structure
To examine whether Scheme IIIA improved the descriptions of the kinetic structure, the 2-D dwell-time distributions and dependency plots predicted by Scheme IIIA were plotted in Fig. 10 for channel B06. During the fitting, the mean lifetimes of the three secondary closed states were made identical by constraining the rate constants for the transitions from the three secondary closed states to the open states to have the same value. Thus, the secondary closed states in Scheme IIIA might be expected to add one additional closed component.
For all five examined channels, Scheme IIIA gave a better description of the Ca2+-dependent kinetics than Scheme III, as indicated in Table I by both the normalized likelihood ratio values, NLR1000, and the rankings, R, by the Schwarz criterion, which applies a penalty for additional free parameters (see Eqs. 3 and 4 in materials and methods).
As would be expected from the improved NLR values and rankings, the plots of the kinetic structure predicted by Scheme IIIA more closely approximated the experimental data, as can be seen by comparing the kinetic structures predicted by Schemes III and IIIA in Figs. 8 and 10, respectively, to the observed and idealized kinetic structures in Fig. 7, B and C. Scheme IIIA predicted the observed excess of long open intervals adjacent to brief closed intervals (#4), which was underpredicted by Scheme III. Scheme IIIA also predicted the magnitude of the deficit of brief open intervals adjacent to brief closed intervals better than Scheme SIII, which overpredicted the deficit (#1). Scheme IIIA also predicted the magnitude of the excess of brief open intervals adjacent to intermediate closed intervals (#2) better than Scheme SIII, although the excess was still overpredicted by Scheme IIIA. Thus, although Scheme IIIA gave a greatly improved description of the kinetic structure over Scheme III (the NLR1000 was increased an average of ∼1,000-fold for the five examined channels), Scheme IIIA still did not capture all of the features.
A measure of the less than perfect description of the gating kinetics by Scheme IIIA is given by the NLR1000. The values of the NLR1000 for Scheme IIIA ranged from 8.9 × 10−7 for channel B14 to 0.0019 for channel M25 (Table I). A value of 1.0 for the NLR1000 indicates that the given kinetic scheme describes the observed distributions as well as theoretically possible for a discrete state Markov model. A value of 0.0019 indicates that the probability that the observed 2-D dwell-time distributions were generated by the given kinetic scheme is only 0.19% per 1,000 interval pairs of the probability that the observed data were derived from the theoretical best description of the distributions. This suggests that Scheme IIIA is still too simple.
Interestingly, an NLR1000 as small as 10−7 can give, by visual inspection, reasonable descriptions of the 2-D dwell-time distributions. The reason for this can be seen by calculating the NLR1 for a single interval pair from the NLR1000 with:
An NLR1000 of 10−7 translates to an NLR1 of 0.984 per interval pair, indicating an average difference in likelihood of 1.6% per interval pair between the predicted and theoretical best likelihood estimates, whereas an NLR1000 of 0.0019 translates to an NLR1 of 0.994 per interval pair. Thus, small differences between the observed and predicted distributions can reduce the NLR1000 far below 1.0. (Note that NLR1 gives average likelihood differences and not average errors, which tend to be greater.)
The small differences in likelihood per interval pair between the predicted and theoretical best likelihood estimates suggests that Scheme IIIA should generate single-channel current records like those from the experimental data. That this is the case is shown in Fig. 11, which presents simulated single-channel currents predicted by Scheme IIIA for channel B06 for comparison to the experimental records in Figs. 7,A and 9. The simulated records include noise and filtering equivalent to that of the experimental data. Due to stochastic variation, the simulated and experimental current records could never look identical, even if the true gating mechanism were known, but the simulated data do capture the kinetic features of the experimental records, including the calcium dependence (Fig. 11,A) and the occasional observation of adjacent brief open and closed intervals (Fig. 11 B).
The rate constants for Schemes III and IIIA are indicated in Fig. 12, A and B, for five examined channels. In general, the values of the specific rate constants were similar for the five channels, although there were notable exceptions. For example, for Schemes III and IIIA, rates k1−4 and k4−1 had values of ∼100 s−1 and 2,000 s−1, respectively, for some of the channels, and values near zero for the other channels.
Since the kinetic structure predicted by Scheme IIIA did not approach the theoretical best description of the Ca2+-dependent gating, a more complex version of Scheme IIIA was considered. Scheme IIIA, with three open and eight closed states, would generate three open and six closed components in the dwell-time distributions because the three secondary closed states have the same lifetimes and would be detected as a single component. Yet, estimates of the minimal number of components for BK channels suggest up to four open and seven to eight closed components (McManus and Magleby, 1988, 1991), suggesting that there may be more detectable states than in Scheme IIIA. Thus, it is possible that the secondary closed states have different lifetimes, since each secondary closed state has a different number of Ca2+ bound, which could alter its stability.
To test this possibility, the lifetimes of the secondary closed states were no longer constrained to be identical during the fitting. The results are indicated under Scheme IIIA′ in Table I, where Scheme IIIA′ improved the NLR for all five channels and ranked higher than Scheme IIIA for four of the five channels. In general, for Scheme IIIA′, the lifetimes of the three secondary closed states differed approximately twofold, but were still brief. For example, the lifetimes of C9, C10, and C11 for Scheme IIIA′ for channel B06 were 0.036, 0.038, and 0.066 ms, respectively.
We also examined whether it was necessary to have three secondary closed states in Scheme IIIA′ or whether two would be sufficient. Deleting the leftmost secondary closed state in Scheme IIIA′ resulted in a lower ranking for four of the five channels (not shown), suggesting that at least three, rather than two, additional closed states may be required. Modifying Scheme IIIA′ to have two secondary closed states gives a model similar to the model of Wu et al. (1995), except that in the model of Wu et al. (1995) the transitions to the secondary closed states are Ca2+ dependent. We examined whether the model of Wu et al. (1995) with Ca2+- dependent transitions to the secondary closed states would improve the descriptions of the kinetic structure. For all five examined channels, this model ranked below Scheme IIIA (not shown). This observation does not necessarily exclude the model of Wu et al. (1995), as they studied BK channels in cochlea hair cells, which may well have different gating properties than BK channels in cultured rat skeletal muscle, and it is also unclear whether they excluded activity in modes other than normal, as we have done.
While Scheme IIIA′ improved the descriptions of the kinetic structure when compared with Scheme IIIA, the increases were small, so that the predicted kinetic structure was still inadequate, appearing similar to that in Fig. 10. Thus, Scheme IIIA′ was still too simple to account entirely for the observed kinetic structure.
Scheme IIIB with Intermediate Closed States Gave Similar Descriptions of the Kinetic Structure as Scheme IIIA′
The above sections explored the effects of adding needed additional brief closed states as secondary states to form Scheme IIIA, in which the secondary closed states would be consistent with channel block or the closing of secondary gates. However, such additional brief closed states might also arise by other mechanisms. Theoretical considerations of gating for a channel with four subunits in which each subunit can bind an agonist and also undergo a conformational change suggests that there may be large numbers of intermediate closed states within the activation pathway (Eigen, 1968; Fersht, 1977; McManus and Magleby, 1991; Cox et al., 1997). It would be difficult, if not impossible, to define the rate constants for such a large number of additional states, but an examination of whether such models might account for gating was made by exploring whether the addition of a few intermediate closed states would improve the description of the kinetic structure. Scheme IIIB extends Scheme III into a minimal model with the three intermediate closed states C9, C10, and C11. The horizontal transition pathways between the three intermediate states have been excluded to reduce the numbers of parameters in the fitting.
Scheme IIIB, which was fitted without placing constraints on the lifetimes of the intermediate closed states, typically gave likelihood values slightly less than Scheme IIIA′. Since Scheme IIIB had the same number of free parameters as Scheme IIIA′, Scheme IIIB typically ranked, on average, just below Scheme IIIA′ (Table I), and predicted a similar kinetic structure (not shown). Thus, Scheme IIIA′ with secondary closed states and Scheme IIIB with intermediate closed states gave similar descriptions of the data. The rate constants for Scheme IIIB are indicated in Fig. 12 C for five examined channels, where it can be seen that the lifetimes of the intermediate states are brief. Similar descriptions of the data were also found (when examined) for the more complex schemes to be presented below when the secondary closed states were replaced with the same number of intermediate closed states. Hoshi et al. (1994) have previously indicated that intermediate closed states can be difficult to distinguish from secondary closed states for steady state data. Thus, in the following sections, only the models with the secondary closed states will be presented in detail, with the implication that the intermediate state models would give similar results.
More Complex Models with Additional States
Four additional classes of models were considered, both with and without secondary closed states. These models explored the effect of adding additional open states in different positions to Scheme III to give four open components, as can often be observed in the experimental data.
Schemes IV and V allow channel opening with one or no Ca2+ bound. Experimental data suggest that such openings may occur (Meera et al., 1996; Cui et al., 1997). Scheme V is the Monod, Wyman, and Changeux (MWC; Monod et al., 1965) model for allosteric proteins. Many of the examined schemes for the BK channel are subsets of this model (McManus and Magleby, 1991; Cox et al., 1997). Scheme VI adds an open and closed state by assuming that the binding of the second Ca2+ can lead to closed states with different properties, depending on whether the second Ca2+ binds to a subunit adjacent to or diagonal from the subunit with the first bound Ca2+. The open states reached from these two different closed states could have different properties, leading to a fourth open component. Scheme VII adds an open and closed state by assuming that the additional open and closed states are reached by a Ca2+-independent conformational change after the binding of four Ca2+. Finally, the effects of adding a secondary closed state to each of the open states for each scheme was examined to obtain Schemes IVA–VIIA.
Based on the average of the rankings for the five channels, Scheme VIIA was the top ranked scheme, followed by Scheme VII, and then Schemes VIA, IIIA′, IIIB, IVA, IIIA, IV, V (MWC), and III. Schemes with secondary closed states ranked above the same scheme without secondary closed states, even though there was a penalty for the additional free parameters associated with the secondary closed states. Replacing the secondary closed states with intermediate closed states (when examined) gave similar descriptions of the data. Thus, these observations are consistent with models in which additional brief closed states, as either secondary or intermediate closed states, contribute to the gating.
The predicted kinetic structures for the schemes ranked above Scheme IIIA in Table I appeared similar to that for Scheme IIIA shown in Fig. 10, so have not been presented. As with Scheme IIIA, the predicted kinetic structures for the higher ranked schemes captured most of the features of the kinetic structure. The observation that Schemes VII and VIIA were the top ranked schemes suggests that the gating may include at least one open and one closed state on the activation pathway that are reached after a Ca2+-independent conformational change. The NLR1000 for Scheme VIIA (Eq. 3) of 0.0067 for channel B06 (Table I) becomes 0.995 for a single interval pair (Eq. 5), suggesting an average likelihood difference per interval pair of only 0.005 between the description of the gating by Scheme VIIA and the theoretical best description for a discrete state Markov model. For all schemes in Table I except Schemes IIIA′ and IIIB, the lifetimes of the additional brief closed states were constrained to be the same for any particular channel. Allowing the lifetimes of the additional brief closed states to vary independently of one another (as was done for Scheme IIIA′) typically improved the descriptions of the kinetic structure and the rankings over the same models in which the brief lifetimes were constrained.
Interestingly, Scheme VIA, the second ranked scheme for four of the five channels, was the top ranked scheme for the atypical channel M09. The predicted kinetic structure for channel M09 for Scheme VIA is shown in Fig. 13, and captures most of the features for this channel (compare with Figs. 2 and 3). The other schemes gave noticeably worse descriptions for channel M09 (not shown). Thus, channel M09 both looks different in kinetic structure and may have some differences in gating mechanism other than differences in rate constants.
In general, for the various examined schemes with secondary closed states, ∼30–40% of the brief closings with durations <110 μs (flickers) arose from transitions to the primary closed states and 60–70% arose from transitions to the secondary closed states. Thus, in terms of the models with secondary closed states, the majority of flickers arose from the secondary closed states. In models with intermediate closed states, essentially all the flickers arose from the intermediate closed states.
Additional Secondary Closed States Can Improve the Rankings
It has been suggested for GABAA channels (Rogers et al., 1994) and for Shaker K+ channels (Schoppa and Sigworth, 1998) that two closed states may be connected directly to the open states beyond the activation pathway. Consequently, we expanded some of the schemes to include two secondary closed states connected to each open state to determine if this improved the description of the kinetic structure for BK channels. For the fitting, all first secondary closed states connected to each open state were constrained to have identical lifetimes and all second secondary closed states connected to each open state were also constrained to have identical lifetimes. In general, such schemes typically ranked above schemes with only a single secondary closed state connected to each open state, even though there were additional penalties for the additional states. Such schemes also gave some minor visual improvements in the kinetic structure. For example, for channel B06, Scheme VIIA had an NLR1000 of 0.0067 with one secondary closed state per open state, compared with 0.0091 after adding an additional secondary closed state per open state. Such schemes have not been presented in detail because the improvement in the kinetic structure was slight for the additional complexity. For experimental data collected over a wider range of experimental conditions, it may well be necessary to include such additional secondary closed states to account for the kinetic structure.
This study employed comparisons of the observed and predicted kinetic structures together with the quantitative ranking of models by simultaneous maximum likelihood fitting of 2-D dwell-time distributions to assess and develop gating mechanisms for BK channels. These approaches extracted correlation information from the single-channel current records, both visually and through likelihood, to place restrictions on the transition pathways and rates for the various kinetic states. In theory, all possible models could be tried blindly and ranked by likelihood to obtain the most likely schemes, but we have found it more efficient to use the kinetic structures to guide the development of the models and assess the ability of the models to describe the single-channel kinetics.
Similarities and Differences in the Kinetic Structure of BK Channels
Similarities and differences in the gating of different BK channels from the same preparation (cultured rat skeletal muscle) were readily apparent in the plots of kinetic structure. While the basic features of the kinetic structure were the same for the different individual BK channels, there were also channel-specific differences (Figs. 2 and 3). The reason for the channel-specific differences in kinetic structure was not established, but alternative splicing of BK channels can give rise to a large family of functionally diverse channels (Atkinson et al., 1991; Adelman et al., 1992; Butler et al., 1993; Lagrutta et al., 1994; Pallanck and Ganetzky, 1994; Tseng-Crank et al., 1994), so it is possible that the individual channels examined in our study may have had different physical structures.
Other factors might also give rise to the differences in kinetic structure. The association of the modulatory beta subunit with the pore forming alpha subunit of BK channels can change the Ca2+ sensitivity and gating kinetics (McManus et al., 1995; Dworetzky et al., 1996; Meera et al., 1996; Nimigean and Magleby, 1998). The extent, if any, to which differential expression of beta (or other possible) subunits might contribute to the differences in kinetic structure is not known for the examined channels. It is also possible that individual channels may have been in different redox states or levels of phosphorylation, which can alter activity (Reinhart et al., 1991; Bielefeldt and Jackson, 1994; DiChara and Reinhart, 1997; Thuringer and Findlay, 1997). Differences in Ca2+ sensitivity and/or kinetics among BK channels from the same preparation have been described previously (McManus and Magleby, 1991; Sansom and Stockand, 1994; Wu et al., 1996). Some of the differences in kinetic structure might also arise from small differences in the experimental conditions for some of the channels, but, as indicated in the results, this is unlikely to be a major contributing factor.
Simple Models Can Approximate many of the Features of the Kinetic Structure of BK Channels
Most of the models that have been examined to describe the gating of BK channels can be related to the MWC model for binding of ligands to a tetrameric protein (McManus and Magleby, 1991; Cox et al., 1997). This model has five closed and five open states, as indicated in Scheme SV. McManus and Magleby (1991) found that Scheme III, a simplified version of Scheme V with two fewer open states, was the simplest model that could describe the basic features of the Ca2+ activation of BK channels (determined with 1-D dwell-time distributions) over an ∼400-fold range of Po, from very low (0.00137) to intermediate (0.53) levels.
As a more critical test of Scheme III, we examined to what extent it could account for the Ca2+ dependence of the detailed kinetic structure. While Scheme SIII captured the basic features, it underpredicted the number of brief closed intervals adjacent to longer open intervals and overpredicted the number of brief open intervals adjacent to intermediate (1 ms) closed intervals (Fig. 8), suggesting that Scheme III was too simple. Based on differences between the observed and predicted kinetic structures, increasingly more complex models with additional states were then examined. Scheme VIIA was the top-ranked scheme in Table I for four of the five examined channels and the second ranked scheme for the atypical channel M09, which had a markedly different kinetic structure from the other channels. Scheme VIA was the top ranked scheme for channel M09.
These observations suggest that the basic underlying gating mechanism is the same for at least four of the five BK channels examined in detail. On this basis, the differences in kinetic structure among these channels would then arise from some differences in the rate constants (Fig. 12). The atypical channel M09 may well have some differences in gating mechanism. Another possibility is that the actual gating mechanism is more complex than any of the examined schemes (see below). If this more complex mechanism were known, then perhaps it would be top ranked for all channels.
Secondary and Intermediate Closed States
A feature that consistently improved the descriptions of the kinetic structure was the incorporation of additional brief closed states connecting directly to the open states. These additional closed states could be outside the activation pathway as secondary closed states, as in the A versions of the various examined schemes, or within the activation pathways as intermediate closed states, as in Scheme IIIB. It was difficult to distinguish between these two types of models. The secondary and intermediate closed states were of brief duration, with typical mean dwell-times of 30–150 μs. The models with the secondary or intermediate closed states could account for additional features of the kinetic structure that were missing with the simpler models (compare Fig. 10 to 7 and 8). For models with secondary closed states, calculations from the most likely rate constants indicated that ∼30–40% of the brief closings with durations <110 μs (flickers) arose from transitions to closed states within the activation pathway and 60–70% arose from transitions to the secondary closed states. For models with intermediate closed states, essentially all of the flickers involved the intermediate closed states.
Consistent with our findings that there are additional closed states, Cox et al. (1997) have found that the basic MWC model can account for the effects of voltage and Ca2+ on Po for mSlo BK channels, but cannot account for the detailed single-channel kinetics. They suggest that additional states are needed. Wu et al. (1995) have added closed states beyond the activation pathway to account for activity in high Ca2+ for BK channels from hair cells in the cochlea. The two added closed states in the model of Wu et al. (1995) are Ca2+ dependent and had intermediate duration lifetimes, compared with the Ca2+-independent briefer lifetimes of the added closed states in our study. Moss et al. (1998) have found evidence for additional brief closed states in dslo BK channels expressed in Xenopus oocytes.
Closed states beyond the activation pathway (secondary closed states) have been proposed for a number of channels including acetylcholine receptors in BC3H-1 cells (Sine and Steinbach, 1987), GABAA receptor channels (Rogers et al., 1994), Shaker potassium channels (Hoshi et al., 1994; Schoppa and Sigworth, 1998) and NMDA (N-methyl-d-aspartate) receptors (Kleckner and Pallotta, 1995). It would be of interest to determine whether models with intermediate, rather than secondary, closed states might also be consistent with the gating of one or more of these channels.
Possible Mechanisms for Intermediate Closed States in BK Channels
The possibility of intermediate closed states follows directly from proposed mechanisms for conformational changes in tetrameric allosteric proteins, in which each subunit can both bind a ligand and undergo one or more conformational changes (Eigen, 1968; Fersht, 1977; McManus and Magleby, 1991; Cox et al., 1997). Transitions from the unliganded channel with all the subunits in a basal conformation to the fully liganded and activated conformations could pass through intermediate states in which the four subunits are in various combinations of conformation with and without bound ligand. Models for tetrameric proteins without intermediate states, such as the MWC model, have been found to be inadequate for hemoglobin, as Ackers et al. (1992) have presented evidence for intermediate states in this tetramer. Thus, unless the intermediate states are too brief or infrequent to detect, it might be expected that some of the intervals during the gating would arise from intermediate states.
If the properties of one or more of the four alpha subunits forming the pore of the channel differ due to alternative splicing (Atkinson et al., 1991; Adelman et al., 1992; Butler et al., 1994), associations with accessory subunits (McManus et al., 1995; Dworetzky et al., 1996; Meera et al., 1996), and/or various levels of phosphorylation and redox (e.g., Reinhart et al., 1991; Bielefeldt and Jackson, 1994; DiChiara and Reinhart, 1997; Thuringer and Findlay, 1997), then the number of potential intermediate states could be many times greater than predicted on the basis of four subunits.
It has been suggested for drk1 and also mutated Shaker channels that intermediate states can give rise to subconductance levels (Chapman et al., 1997; Zheng and Sigworth, 1997). Hence, if the brief closings (flickers) for BK channels arise from transitions to intermediate states, then it might be expected that brief closings would be to subconductance levels. Consistent with this possibility, the findings of Ferguson et al. (1993) suggest that brief closings (flickers) for BK channels are, on average, partial closings to ∼90–95% of the fully closed conductance level.
If the intermediate states are partially conducting, then, since transitions to and from the open states in Scheme IIIB must pass through the intermediate states, it would be expected that the transitions between the fully open and fully closed current levels in the single-channel records would pass through brief lifetime subconductance levels. This hypothesis is consistent with the observations of Ferguson et al. (1993) that openings from longer duration closed intervals first passed through, on average, a brief (50 μs) subconductance level at the 90–95% closed level before opening fully, and that closings to the longer duration closed intervals first closed, on average, to the 90–95% closed level before closing completely over the next 50 μs.
Possible Mechanisms for Secondary Closed States in BK Channels
While there is a theoretical underpinning for intermediate states within the activation pathway (see above), the basis for potential secondary closed states, those beyond the activation pathway, is less clear, but some possibilities include either blocking ions or secondary gates that obstruct the channel for brief periods of time when the primary gate is open. In either case, if secondary states give rise to the brief closings (flickers), then transitions to the secondary closed states, on average, would need to be associated with subconductance levels to be consistent with the observations of Ferguson et al. (1993) discussed above. (Furthermore, to account for the subconductance levels that can be observed, on average, for transitions between the fully open and fully closed current levels in the single-channel records, some states in the activation pathway would still need to be partially conducting.)
Although some blocking ions such as Ba2+ decrease the conductance to the zero current level without measurable delay (Ferguson et al., 1993), blocking ions can also give partial conductances (Schild and Moczydlowski, 1994; Premkumar and Auerbach, 1996; Shneggenburger and Ascher, 1997). Thus, blocks could serve as a mechanism for secondary closed states with either full or partial closures. While Na+ produces a fast flickery block in BK channels (Yellen, 1984), our experimental solutions contained negligible Na+, so Na+ block can be excluded for secondary closed states.
It is also unlikely that the large numbers of brief closed intervals in our study arose from fast Ca2+ blocks since brief closed intervals do not become more frequent with increasing Ca2+i (Rothberg et al., 1996), and millimolar concentrations of Ca2+i are required to reduce the conductance by apparent screening effects (Ferguson, 1991). Intracellular Ba2+ applied to BK channels can produce long lasting blocks with durations of seconds (Vergara and Latorre, 1983; Rothberg et al., 1996) as well as fast flickery blocks (Sohma et al., 1996; Bello and Magleby, 1998). While an occasional flicker in our records may arise from a Ba2+ block, it is unlikely that an appreciable number do, as the expected submicromolar concentrations of contaminant Ba2+ in our solutions would be three to four orders of magnitude less than the concentration required for appreciable fast Ba2+ blocks, and the mean durations of slow Ba2+ blocks are five orders of magnitude longer than the flickers (Rothberg et al., 1996). Thus, the flickers in our study do not appear to arise from ion block of the channel by Ba2+, Ca2+, or Na+.
Therefore, the flickers might arise from the closing of a secondary gate intrinsic to the channel or the closing of the primary gate to a secondary conformation. Mienville and Clay (1996) have proposed a secondary gate that acts independently of the main gate to describe the flickering in BK channels in embryonic rat telencephalon. They found that decreasing intracellular K+ increased the flicker rate, leading to a decrease in mean open time, suggesting that the open state may be destabilized when the pore is unoccupied by the permeating ion (Swenson and Armstrong, 1981; Demo and Yellen, 1992). Any secondary gates associated with the secondary closed states would have properties that differ from the much slower inactivation gate observed by Solaro and Lingle (1992) for BK channels in rat adrenal chromaffin cells.
Minimal Working Hypothesis for the Gating in the Normal Mode
Of the examined models, Schemes IIIA, IIIA′, and IIIB can be considered as minimal models that capture the major features of the Ca2+ dependence of the gating, as revealed by plots of the kinetic structure. Fig. 14 summarizes the features of the gating for Scheme IIIA with the rate constants in Fig. 12,B. Fig. 14,A plots the rate constants for binding of the first through fourth Ca2+ to the channel. The binding of the second and third Ca2+ to the closed states was highly cooperative, with the second Ca2+ binding rate constant, on average, 9× faster than the first, and the third Ca2+ binding rate constant, on average, ∼100× faster than the first. The fourth Ca2+ binding rate constant was not as fast as the third, but was still ∼80× faster than the first. For a channel with four identical subunits with one Ca2+-binding site per subunit, the Ca2+-binding rate constants per binding site would be 1/4, 1/3, 1/2, and 1/1 of the indicated rate constants in Fig. 14,A, increasing the forward binding rate for the successive bindings of the second, third, and fourth Ca2+ even further. The cooperativity in the binding of successive Ca2+ to the closed states for Schemes IIIA′, IIIB, IV, IVA, V, VA, VII, and VIIA was similar to that shown in Fig. 14 A for Scheme IIIA. Cooperativity requires that the subunits interact to facilitate the binding of additional Ca2+. If the channel opened after binding three Ca2+, the rate constant for binding the fourth Ca2+ to the open state was, depending on the model, typically similar to or several times faster than the rate constant for the binding of the fourth Ca2+ to the closed state.
Fig. 14, B and C, plots the mean lifetimes of the open and closed states for Scheme IIIA. In general, increasing the number of bound Ca2+ increased the stability (lifetimes) of the open states and decreased the lifetimes of the longer closed states. This was also generally the case for the other schemes, including Scheme IIIB. It is this general inverse relationship between the durations of the longer closed states (including compound closed states) and the open states (including compound open states) that gives rise to the characteristic saddle shape of the dependency plots.
In Scheme IIIA, the main functional gating pathway was typically C8-C7-C6-C5-O2-O1, with brief transitions to the secondary closed states C9 and C10 interrupting the openings, but having little effect on Po. For Scheme VIIA, the main gating pathway was typically C10-C9-C8-C7-O3-O2-O1, with brief transitions to the secondary closed states C11, C12, and C13 interrupting the openings. The equilibrium occupancies of the various states as a function of Ca2+i are plotted in Fig. 14, D and E. Ca2+i drives the channel towards the closed states C4 and C9 and the open state O1.
The Ca2+-binding rates, state lifetimes, and equilibrium occupancies for Scheme IIIB were similar to those shown in Fig. 14 for IIIA, except that the equilibrium occupancies of the intermediate closed states C9, C10, and C11 were typically about twice those for the secondary closed states C9, C10, and C11 for Scheme IIIA. For Scheme IIIB, the main functional gating pathway was typically C8-C7-C6-C5-C10-O2-O1, with brief transitions to the intermediate closed states C9, C10, and C11 interrupting the openings, but having little effect on Po.
Interestingly, while some of the models we considered are similar to those considered by Cox et al. (1997) for the gating of mslo, there is little agreement between the two studies for the estimated parameters in the considered models. This lack of agreement is not surprising. Cox et al. (1997) fit macroscopic currents over wide ranges of Ca2+i and voltage, but did not include single-channel kinetic data, and consequently could not account for the single-channel kinetics. We fit single-channel data over a lesser range of Ca2+i at a single voltage (+30 mV), and consequently could account for single channel kinetics, but not voltage. The models from both studies summarize the important observations for the examined experimental conditions, while placing some restrictions on mechanism. These models can serve as minimal models for future studies.
The Top Ranked Models Are Still Too Simple to Account for Normal Activity
BK channels can open at nominally zero Ca2+i (Rothberg and Magleby, 1996; Cox et al., 1997; Meera et al., 1997). Of those considered, only Schemes V and VA would open in the absence of Ca2+i. Thus, to account for the gating at very low Ca2+i, the more complex Scheme VA, or a version of this scheme in which the secondary closed states are replaced with intermediate closed states, could serve as a minimal model. Open states could also be added to the left-most closed states in the other schemes to allow opening in the absence of Ca2+i.
Alternatively, while the horizontal steps in our considered schemes are listed as Ca2+-binding steps, each step may actually represent a Ca2+-facilitated conformational change. If the conformational changes can occur spontaneously, then this could contribute to some of the gating at very low Ca2+i. It is also possible that the various cations in the bathing solution might facilitate spontaneous conformational changes, but at a much lower efficacy than Ca2+. If this is the case, then a spontaneous transition rate for the forward directed horizontal transitions could be added to the models.
While the top ranked models gave excellent descriptions of the kinetic structure over the examined range of Ca2+i, they still did not describe the data as well as theoretically possible for a discrete state Markov model. Although the differences between the predicted and observed kinetics were small, such differences are important, as they may give additional insight into the mechanism. Small amounts of drift in the single channel data, minor errors in the estimated Ca2+ concentrations, and errors associated with the assumptions used to correct for missed events (Crouzy and Sigworth, 1990) could contribute some error. The assumptions of discrete state models with rate constants that do not change in time may also contribute some error. Although the gating of BK channels appears consistent with these simple assumptions (McManus and Magleby, 1989; Petracchi et al., 1991), conformational changes in proteins would be expected to be more complex (see discussion in McManus and Magleby, 1989).
Alternatively, the predictions of the models may be less than complete because the examined models are still too simple, having too few kinetic states. As mentioned above, theoretical considerations predict far more complex models for BK channels than the minimal models examined in this study (McManus and Magleby, 1991; Cox et al., 1997). If each subunit has two Ca2+ binding sites (Golowasch et al., 1986; Schreiber and Salkoff, 1997), then the complexity of the models would increase even further.
The trouble with more complicated models is that it becomes increasingly more difficult to define the parameters as more states are added. Consequently, even though the more complex models may give significantly better descriptions of the data, the number of parameters that are not well defined typically increases with complexity, and the rankings of more complex models can become inconsistent among channels with so many free parameters, so that it may not be clear which of the more complex models should be used.
For example, we found that models with two secondary closed states connected to each open state typically improved the description of the kinetic structures and also ranked higher than models with only one secondary closed state connected to each open state. Based on this observation, it would also be expected that models with additional intermediate closed states, or with both intermediate and secondary closed states, would rank higher than models with only secondary or intermediate closed states alone. However, rankings of the more complex models that were examined were no longer as consistent among the different channels as for the simpler models. Another complicating factor is that more complex schemes, even if correct, may be ranked lower by the Schwarz (1978) criteria if the additional states do not significantly improve the fit. Thus, the detailed testing of more complex models seems best delayed until there are more data to constrain the models, such as detailed single-channel currents obtained over wider ranges of Ca2+i than used in these experiments and also over a range of voltages as well as following step changes in voltage. Such data should provide further constraints on the parameters.
In conclusion, our study describes and presents minimal models to account for the kinetic structure of BK channels. These models have additional brief closed states, either in the activation pathway as intermediate states or beyond the activation pathway as secondary states. It is these additional brief closed states that give rise to the majority of the brief closed intervals (flickers) in the gating.
This work was supported in part by grants from the National Institutes of Health (AR-32805, NS-30584, NS-007044) and the Muscular Dystrophy Association.
Abbreviations used in this paper
Address correspondence to Dr. Karl L. Magleby or Dr. Brad S. Rothberg, Department of Physiology and Biophysics, University of Miami School of Medicine, P.O. Box 016430, R-430, Miami, FL 33101-6430. Fax: 305-243-6898; E-mail: firstname.lastname@example.org