| HOME | HELP | FEEDBACK | SUBSCRIPTIONS | ARCHIVE | SEARCH | TABLE OF CONTENTS |

* Department of Physiology and Pharmacology, State University of New York Downstate Medical Center, Brooklyn, New York 11203; and
University of Wales, Singleton Park, Swansea SA2 8PP, United Kingdom
Correspondence: Address reprint requests to James J. Celentano, E-mail: jimc10{at}aol.com.
| ABSTRACT |
|---|
|
|
|---|
1ß2
2S, 1 mM GABA, >2.5 s), a model with two open states and three desensitized states is favored. When applied to simulated data generated using a consensus model, CVF leads to reasonable parameter estimates and accurately discriminates between this and other models. | INTRODUCTION |
|---|
|
|
|---|
To date, most studies have focused on measuring surrogate parameters such as relaxation time constants or mean sojourn times. From this work, several general features of GABAR function have become apparent. Open and closed time distributions of single-channel events under steady-state conditions are best fit by multiple components (Macdonald et al., 1989
; Weiss and Magleby, 1989
). From this, it is clear that the receptor is capable of visiting multiple distinct open and closed states. Because these distributions are derived from binned data, however, a great deal of information about the time correlations is lost. Studies examining the correlations between adjacent events (Ball and Sansom, 1989
; Horn and Lange, 1983
; Magleby and Weiss, 1990
; Qin et al., 1996
, 1997
; Sakman and Neher, 1995
; Weiss and Magleby, 1989
) as well as the hidden Markov (Qin et al., 2000a
,b
) and the HJCFIT (Colquhoun et al., 2003
) methods begin to take advantage of this type of information. The most comprehensive kinetic model, based on single-channel and macroscopic GABAR data, has a total of three open states, nine closed states, and three desensitized states (Haas and Macdonald, 1999
). Two fundamental questions remain at least partially unanswered: how many of theses states are relevant to the description of synaptic function and what is the exact arrangement (i.e., connectivity) of the relevant states?
Perturbation studies using agonist concentration jumps examine the response of a population of GABARs to large and rapid changes in agonist concentration. This experimental paradigm more closely mimics synaptic events and the resulting traces are typically fit to exponential functions by sum-of-squares minimization. This method also ignores time correlations (see Appendix I) and generates surrogate values (time constants and magnitudes) whose relationship to the underlying molecular mechanisms is unclear. The recent application of first latency analysis to perturbation results begins to examine the arrangement of states and suggests that desensitization can precede channel opening (Burkat et al., 2001
).
What is needed is a method that can fit macroscopic perturbation data directly to a particular kinetic model by taking into account both the magnitude of the current at each time point as well as the correlations among the different points. The method described here, referred to as covariance fitting (CVF), accomplishes this by using the Q-matrix technique to calculate the correlations predicted by various kinetic models. Coupled with the general multivariate likelihood function and likelihood maximization using a variable metric algorithm, this method accurately estimates kinetic parameters and reliably discriminates between models with up to 11 free parameters.
| METHODS |
|---|
|
|
|---|
For models with a single open state: 1D, 2D, 3D, 4D indicate 1, 2, 3, or 4, desensitized states, respectively. The number of desensitized states on each side of the gating isomerization has a powerful influence on CVF results. To illustrate this feature of each model, a Roman numeral is used to indicate the number of desensitized states on the closed side of the gating isomerization. Thus, model 3D-I has a total of three desensitized states with one connected to the closed state and two connected to the open state. The arrangement of two or more desensitized states on the same side of the gating isomerization can be either linear (ln) or branched (br). Unless otherwise stated, all such models are branched, meaning all desensitized states are connected directly to either the open or closed state. 1D-C refers to a cyclic model in which the D state is connected to both the closed and open states. 3D-C refers to a cyclic variation of 3D-I in which a transition is allowed between the D state connected to the closed state and one of the D states connected to the open state. This nomenclature is not meant to be comprehensive, but uniquely describes the single open-state models studied here.
No systematic approach is used for naming the limited number of two and three open-state models studied here. Unless otherwise indicated, all multiopen-state models have three desensitized states. See legend to Fig. 7 for the abbreviations used to describe multiopen-state models.
|
) and beta (ß) refer, respectively, to leaving or entering an open state from either the closed state or from a more proximal open state.
Generation of simulated data
Simulated stochastic data are generated using various kinetic models. The agonist binding step has been omitted and simulations are done as if GABA association is instantaneous. At time zero, all receptors start in the closed liganded state (C). Two different methods are used to choose subsequent dwell times. The bulk of the simulated data sets are generated using the QUB software (Research Foundation, State University of New York, Buffalo, NY) (Qin et al., 1996
, 1997
) that determines dwell times directly from the results of a random number generator with an exponential distribution and an average equal to the inverse of the exiting rate constant. A limited number of simulated data sets are generated using the results of a random number generator with a uniform distribution. Here, a random number between 0 and 1 is generated for each 0.01-ms time interval. The dwell time ends and a transition occurs when the random number falls within a critical region. The size of the critical region is proportional to the probability that a transition will occur during a 0.01-ms period. The fastest rate constant used (3000 s1) has a transition probability of only 0.03 during any given interval. No differences are observed between the two methods in any subsequent analysis.
Smooth simulations are done either by numerically solving the differential equations using the Runge-Kutta method, or by using the Q-matrix technique (Eq. 3, Appendix I). The values generated by the two different methods for the same model differ by <0.3%.
Zero-time determination
Because the concentration jump in an actual experiment is not instantaneous, the most appropriate time point to be used as time 0 must be determined empirically. This is done by fitting the rising phase (beginning when the current reaches half the maximum) and the initial decay phase to the following formula using sum-of-squares minimization (Prizm, GraphPad, San Diego, CA):
![]() |
r is the time constant of the rising phase,
d is the time constant of the fast phase of desensitization, Id is the magnitude of the fast phase of desensitization, Iss is the magnitude of the subsequent phases of desensitization plus the steady state, and
t is adjusted as a free parameter. For subsequent analysis, each time point is adjusted by the value of
t arrived at using the above fit.
Data point selection
Because CVF requires generation of a covariance matrix (see Appendix I), it is typically impossible to analyze an entire data set. A 3000-point data set, for example, would require a matrix with 9,000,000 elements. Data sets of 100200 points are, therefore, selected from the raw sets. Selecting points evenly spaced in time (every 30 ms for a 3-s trace) would sacrifice a great deal of information during the steep initial phase of desensitization. Here, points are selected so that there is a roughly even distribution along the idealized path of desensitization (Fig. 1).
|
|
Inorm = I0 If), where I0 is the current extrapolated to t = 0, and If is the current predicted by the exponential equation at the final time point.
Inorm is used to normalize each current value. Time is normalized by dividing each value by
tnorm x tsf, where
tnorm is the final time point and tsf is a timescale factor that can be adjusted to favor earlier or later points. Unless stated otherwise, selection is done evenly along the path of desensitization by setting tsf = 1. When tsf > 1, earlier time points are favored and when tsf < 1, later time points are favored.
The normalized distance along the idealized path of desensitization is calculated in segments (Si) using the Pythagorean theorem (see Fig. 1).
Ii is determined by evaluating the exponential equation for each time point associated with the raw data and
ti is the interval between each time point. As the segments are added, each time point is associated with a distance along the path of desensitization. The total distance is calculated by summing the normalized segments. The total distance along the idealized path is divided by the number of points desired (100200) minus 1, giving a step distance. Finally, the selected data set is created using the relationship between the time and the distance along the path.
The first point in the raw data set with a current greater than one-half the peak current is included as the first point in the selected data set. Each subsequent point included is the first point past the next step distance. If a step distance is not associated with a point, the next point is included and the step distance recalculated based on the remaining distance and the number of points yet to be included.
Parameter estimation
CVF estimates kinetic parameters using the principals of maximum likelihood (ML) (see Appendix I) and a variable metric algorithm. A clear description of the variable metric algorithm can be found in Acton (1970)
, however the formula for calculating the Hessian matrix differs from that shown in Press et al. (1992)
. This study uses the formula presented in Press et al. (1992)
. First and second derivatives of the likelihood with respect to each parameter are calculated numerically. Maximization along a chosen trajectory is done by parabolic interpolation (Press et al., 1992
). Convergence is reached when subsequent iterations change the parameters by a factor of <0.0005. This allows CVF to routinely identify differences in L[LH] as small as 0.001 (see Fig. 5 B). A detailed description is given in Appendix II of the different methods used to generate initial values for CVF. When performed using a Pentium IV 2.4 GHz processor (Intel, Santa Clara, CA), fitting a single data set with 200 points to a model with 11 free parameters typically requires between 12 and 36 h.
|
Model discrimination
In addition to identifying parameter estimates, the ML method generates a maximum likelihood estimate (MLE) that is used to discriminate between different kinetic models. Two models are compared by fitting the same data set to each model and calculating the ratio of their respective MLEs, the likelihood ratio (LR). As a test of statistical significance, it is well established that no test has a higher power than the LR test (Mood et al., 1974
; Rao, 1973
). To apply the LR test, it is necessary to know the expected distribution of the natural logarithm of the LR (L[LR]) when one of the two models (the null model) is correct. Highly improbable values of L[LR] lead to rejection of the null model.
When comparing models with the same number of free parameters, the model with the higher MLE is generally favored. When comparing models with a different number of free parameters, the more complex model often has a higher MLE simply by virtue of its greater flexibility. To legitimately reject the simple model, it is necessary to show that the higher MLE is significantly greater than that predicted by random chance while taking into account the extra free parameters. This requires knowing the distribution of the L[LR] when the simple model is correct.
If two models are identical except for the existence of an additional state, the simple model is said to be a submodel of the more complex model. The two are identical when the value of the forward rate constant is 0 or the value of the reverse rate constant is infinite (i.e., the receptor never enters or remains in the additional state). When comparing a complex model and a submodel, it is well known that, under specific conditions, 2 x L[LR] has a
2-distribution with the number of degrees of freedom (df) equal to the number of unnecessary free parameters (Bishop et al., 1975
; Horn, 1987
; Knight, 2000
; Shao, 1999
). This relationship holds when the submodel is the correct model, the parameter estimates are approximately normally distributed, and there are no constraints on the free parameters in the neighborhood of their correct values (Chernoff, 1954
; Feder, 1968
; Kudo and Arbor, 1963
; Nuesch, 1965
; Perlman, 1969
; Self and Liang, 1987
; Shapiro, 1985
). Kinetic parameters, however, can only be
0, which leads to a violation of this last condition. The correct value for the parameter associated with entry into the unnecessary state is 0 when the simple model is correct. As a result, 2 x L[LR] when comparing a simple model to one containing an additional unnecessary state will not necessarily follow a
2-distribution. In this study, the results of fitting simulated data sets to different kinetic models is used to empirically describe the L[LR] distribution when a correct model is compared to various alternate models. Because L[LR]s are not normally distributed and different comparisons have different L[LR] distributions, tables reporting L[LR]s give the median values as well as the average and range of the values.
Statistics and p-values
Pairs of kinetic models (A and B) are compared using the LR test. Model A is rejected in favor of B if it can be shown that the probability of the observed L[LR] distribution is <0.05 when model A is known to be correct. In cases where A and B have the same number of states and free parameters, this will be done using a simple binomial distribution. If the null hypothesis is that the two models cannot be distinguished, meaning the L[LR] is distributed evenly around 0, the probability of observing 10 out of 10 sets favoring model B, for example, would be (1/2)10, corresponding to p < 0.001. The probability of observing <10 out of 10 sets is calculated using the general binomial equation (Bevington and Robinson, 1992
).
Calculation of p-values using a binomial distribution is nonparametric and does not take into account the actual magnitude of the L[LR]s. In some cases where A and B have the same number of free parameters, the p-value does not favor rejection of one model over the other. In these cases the sum of the natural logarithm of the MLEs (L[LH]s) can sometimes be used to decide which model is favored (see below). No attempt is made here to determine a p-value based on the magnitude of the sum of the L[LH]s. Models supported by the L[LR] distribution and a p < 0.05 are referred to as "heavily" favored, and models supported by the sum of the L[LH]s alone are referred to as "slightly" favored.
Simulated data are analyzed to determine the effect on L[LR] of the inclusion of an unnecessary desensitized state. Results indicate that when comparing the correct model to one with an unnecessary state, the median L[LR] is
0.1 meaning
50% of the L[LR]s are between 0.0 and 0.1 and 50% are >0.1. Simple models can therefore be rejected based on a p-value calculated using a simple binomial distribution. Unless otherwise stated, all p-values reported here are calculated using the binomial equation. For comparisons between complex and submodels, the null hypothesis is that the probability of L[LR] being >0.1 is 50%. When >50% of the L[LR]s are <0.1, the p-value is reported as greater than the value when exactly 50% are <0.1. For comparisons between models with the same number of states, the null hypothesis is that the median L[LR] is 0 meaning the probability of L[LR] being >0 is 50%. When the median falls within 0.001 of 0, the p-value is reported as greater than the value when exactly 50% are <0.
Unless stated otherwise, all averages are expressed mean ± SE (
/N1/2). All logarithms are natural logarithms.
Recombinant receptors
HEK 293 cells are grown in MEM supplemented with 10% fetal bovine serum. Transfection with DNA coding for the
1ß2
2s subunits of the GABAA receptor (1:1:1.2) is done by CaPO4 precipitation (Chen and Okayama, 1987
). To facilitate the identification of transfected cells, DNA coding for the CD8 antigen is included. Transfected cells are identified using conventional light microscopy as those that bind anti-CD8 coated beads (Jurman et al., 1994
).
Electrophysiology
Recordings are done at room temperature (2022°C) in standard recording buffer (in mM, 150 NaCl, 4 KCl, 1 MgCl2, 1 CaCl2, 10 HEPES, pH 7.2 with NaOH) on the stage of an inverted phase contrast microscope enclosed in a Faraday cage. Patch electrodes (813 M
) pulled from thin-walled borosilicate glass are filled with electrode buffer (in mM, 3 NaCl, 140 KCl, 11 EGTA, 10 HEPES, 4 MgATP, pH 7.2 with KOH).
GABA (1 mM) is applied to voltage-clamped (from 50 to 70 mV) outside-out patches (2.53 s) using a
-tube driven by a piezoelectric translator. The 1090% exchange rate is <1 ms. Currents are recorded using a Dagan 3900A integrating patch-clamp amplifier (Minneapolis, MN), low-pass filtered (2 kHz) and digitized at 0.25 ms per point, and stored on disk for further analysis offline. Because the rise time of the filter (0.17 ms) is shorter that the sampling interval, no correction of the covariance matrix is made for the effects of the filter.
| RESULTS |
|---|
|
|
|---|
The LR test is used to make three general types of model comparisons: models with a different number of states (or free parameters), models with the same number but a different arrangement of states, and models with the same number and arrangement of states but with parameter values corresponding to different patterns of desensitization (see below). Unless otherwise stated, the following conventions are followed. When comparing models with a different number of states, the LR will be the MLE associated with the complex model divided by the MLE associated with the simple model. When comparing models with the same number of states but with a different arrangement of states or a different pattern of desensitization, the LR will be the MLE associated with the correct model divided by the MLE associated with the alternate model. Scatter plots showing the L[LR] distributions are labeled "Number", "Arrangement", or "Pattern" to indicate the type of comparison being made.
Simulated 1D data
The simple model 1D-I is used to generate a total of 80 data sets of 100 channels each 300 ms in duration and 0.1 ms per point (see Fig. 1 for a sample plot; see Table 1 for simulation parameters). The kinetic parameters are chosen to give a desensitization time constant of
30 ms.
The first 40 data sets were originally generated to determine the distribution of exponential parameters derived from fitting by sum-of-squares minimization. When the raw data were fit to single and double exponential equations, however, the F-test incorrectly favored the latter equation for 11 of the 40 sets. Furthermore, the p-value for 8 of the 11 is <0.0001. An accurate method for discriminating between the two equations should favor the incorrect equation, on average, in only two of the 40 sets (p-value <0.05). The F-test and sum-of-squares minimization assume the variance is uniform and uncorrelated. CVF was developed both to overcome the shortfalls of sum-of-squares minimization and to take advantage of the information contained in the covariance of ion channel data.
For each of the 80 simulated data sets, 100 points are selected evenly along the idealized path of desensitization (Fig. 1) and fit using CVF to the correct model (1D-I). The resulting kinetic parameters are approximately normally distributed (Fig. 2 A) as predicted by ML theory. For three of four parameters, the parameter estimate is within two standard errors (SEs) of the true value and all four estimates are within 10% of the true value (Table 1). The L[LH]s are also approximately normally distributed, as predicted by ML theory (Fig. 2 B). When the same 80 data sets are fit to the same model using wSS (see Methods), the SEs are larger and only two of the four estimates are within 10% of the true value (Table 1).
|
When each of the 80 data sets generated using 1D-I are fit to 2D-II, the L[LH]s are approximately normally distributed (Fig. 2 B). The L[LR] associated with comparing 2D-II and 1D-I is determined for each set. Thirty-seven sets have a L[LR] of <0.1, including 25 with a L[LR] of <0.001. The distribution of 2 x L[LR] values differs substantially from that of a
2-distribution with 2 df (degrees of freedom) (Fig. 2 B). Fig. 3 A shows a scatter plot of the individual determinations of L[LR] and the values are summarized in Table 2. The average L[LR] value is 0.546, which differs from 1, which is the value expected if 2 x L[LR] has a
2-distribution with 2 df.
|
|
|
|
To test the ability of CVF to discriminate between models with the same number but a different arrangement of states, 30 of the 80 sets generated using 1D-I are fit to a model with the desensitized state connected to the open state (1D-0). All but three have positive L[LR]s (Fig. 3 B). Compared to a null hypothesis in which CVF has no discriminatory power and L[LR] is distributed evenly around 0, 27 out of 30 values above 0 are associated with a p-value of <105, based on the binomial distribution (see Methods). The values illustrated in Fig. 3 B are summarized in Table 3.
Ten data sets are generated using model 1D-0 and fit to both models. Nine of 10 sets had positive L[LR]s, favoring the correct model (Fig. 3 B, Table 3). Compared to a null hypothesis in which L[LR] is distributed around zero, this is associated with a p-value of <0.01. Despite having equal numbers of free parameters, models 1D-I and 1D-0 result in data sets that can be reliably distinguished using CVF and the LR test. In both cases, the correct model is heavily favored over the incorrect model.
When 1D-I and 1D-0 data sets are fit to both models using wSS, 1D-0 is favored in both cases (Table 3), suggesting that the ability of CVF to accurately discriminate between models with a different arrangement of states comes from its ability to extract correlation information from the data.
Simulated 2D data
When simulated data sets generated using four different 2D models (10 sets each) are analyzed by CVF, similar results are observed as with 1D data sets. Twenty-one of 24 parameter estimates are within two SEs of the correct value (Table 4). For all four models, 10 of 10 L[LR]s are >0.1 when fitting to 1D-I and well over half of the L[LR]s are <0.1 when fitting to 3D models (Fig. 4 A, Table 2).
|
Because the two desensitized states are indistinguishable, data generated using model 2D-I (D1==C==O==D2) fit to the correct model with at least two distinct local maxima depending on which of the two desensitized states is responsible for each of the two phases of desensitization. To define their positions in relation to the gating isomerization, desensitized states connected to the closed state (such as D1 above) will be referred to as "proximal" states (p) and states connected to the open state (such as D2 above) will be referred to as "distal" states (d). 2D-I models can then be divided into 2D-I,pd and 2D-I,dp to indicate the pattern of desensitization. The former term refers to models in which the proximal state is responsible for the fast phase and the latter refers to models in which the distal state is responsible for the fast phase. CVF can identify the local maximum associated with the correct pattern of desensitization, when the rates of the two phases of desensitization are sufficiently different (Fig. 4 C, Table 5). As above, only seven of the 10 L[LR]s for 2D-I,3-700,dp are clearly above 0 when comparing the dp and pd patterns. The results are also said to slightly favor the correct model based on the magnitude rather than the distribution of the L[LR]s (p = 0.117). The total of all 10 L[LH]s for the correct model exceeds that for the pd pattern by 24.8 log units.
|
|
To confirm that CVF can identify four desensitized states, the model 4D-II is used to generate a series of data sets with 100 channels each (see Table 6 for parameters). Initial results from data sets of 100 evenly selected points give very small values for L[LR] when comparing 4D-II and 3D-I (not shown). A series of eight data sets are then generated with 200 channels each and 200 points are evenly selected. When these data are fit to 4D-II and 3D-I, only one set had a L[LR] of <0.1 (Table 2). The magnitudes of the L[LR]s are, however, not as great as those observed when 3D data are fit to 3D and 2D models.
To test the ability of CVF to discriminate between different arrangements of states, one of the two groups of 3D-I data sets is fit to a series of 3D models (Fig. 5, Table 3). When compared to the correct model, 3D-II is slightly disfavored (based on total L[LH] alone), whereas both of the entirely asymmetric models (3D-III and 3D-0) are heavily disfavored (p < 0.01 and p < 0.001, respectively). Fitting using wSS fails to discriminate between 3D-I and 3D-0 (Table 3).
In fitting to either 3D-I or 3D-II, it is impossible to discriminate between branched and linear forms. The value of the L[LR]s for these comparisons are rarely >0.001 above or below 0. This is surprising given that, in the case of 3D-I, the comparison involves models with differing numbers of gateway states (states connected directly to an open state). The number of gateway states dramatically affects the correlations observed when fitting adjacent events recorded from single-channel data (Ball et al., 1989
; Colquhoun and Hawkes, 1987
). Fig. 5 B illustrates the distribution of L[LR]s when it is impossible to discriminate between models, when one model is slightly favored, and when one model is heavily favored over another. This figure also shows that CVF can regularly identify differences in L[LH] as small as 0.001.
The 3D-I models used to generate the two groups of simulated data sets differ from one another in the pattern of desensitization. One has parameters that give a pattern of desensitization in which the proximal state is responsible for the fast phase and the two distal states are responsible for the intermediate and slow phases (pdd). The other group has parameters that give a pattern of desensitization in which the distal states are responsible for the fast and intermediate phases and the proximal state is responsible for the slow phase (ddp). CVF heavily favors the correct pattern for both 3D-I,pdd and 3D-I,ddp data sets (Table 5). When applied to 3D-I,pdd, wSS fails to discriminate between the two patterns of desensitization (Table 5). Taken together, these results demonstrate that CVF can reliably extract kinetic parameters and discriminate between a variety of kinetic models with a single open state and up to three desensitized states.
Model discrimination with experimental data
GABAR data are recorded from six different outside-out patches as described in Methods (see Fig. 6 for sample trace) and analysis proceeds in several phases. The standard deviation of the baseline current is squared and used as an estimate of the background variance. For CVF, this value is added to the diagonal elements of the covariance matrix. The baseline current noise has an average covariance of 0.45 pA2 for an interval of 0.25 ms. No effort is made to correct for this degree of correlation because the predicted covariance from channel activity for the same time interval is 21 pA2 based on the final consensus model (see below). No effort is made to correct for open-channel noise, which is expected to have a standard deviation of
3% of the single-channel current (Sigworth, 1985
). The standard deviation of the background current has an average of
56% of a single-channel current.
|
1520%, with a compensatory increase in the former. For this reason, the conductance per channel is fixed at 28 pS in the analysis of GABAR data.
Optimization of point selection
The results of each exponential fit are also used in selecting points along the idealized path of desensitization, as described in Methods, and to generate initial starting values for CVF (see Appendix II). Initial CVF results are used to create consensus models that are then used to generate simulated data. Analysis of the simulated data leads to refinement of the point selection process, CVF is reapplied to the six reselected GABAR data sets and the entire process repeated. The intermediate results from the first two selected GABAR data sets are described briefly, followed by a detailed description of the final analysis using the optimally selected points.
The process begins by selecting from each of the six raw traces 100 points distributed evenly along the path of desensitization and fitting each set to 11 different single open-state models (as shown for 3D-I data in Fig. 5). The number of channels is constrained to the estimate made as described above. Results on the number and arrangement of states are similar to those observed when 3D-I,pdd data are analyzed. 3D models are heavily favored over 2D models whereas at least half the 4D models had L[LR]s <0.1. 3D-I is slightly favored over 3D-II whereas the fully asymmetric models, 3D-0 and 3D-III, are heavily disfavored. For both 3D-I and 3D-II, most comparisons between branched and linear forms have L[LR]s within 0.001 of zero. The ddp pattern of desensitization is favored over pdd.
A key feature of the ddp pattern of desensitization is a high value for the forward rate constant associated with the proximal desensitized state and a low value for the corresponding reverse rate constant (see Table 6). This predicts that a fraction of the receptors desensitize by entering the proximal state without ever opening. Entry of the remaining receptors into the proximal desensitized state, which is thermodynamically favored, is slow because so few receptors (<5%) are in the closed state once the first gating isomerization has reached equilibrium (<1 ms). The fast and intermediate phases of desensitization reflect equilibration with the two distal states.
The estimate of channel number based on fitting to exponential equations may be unreliable if a significant proportion of the channels desensitize without opening, as predicted by the ddp pattern of desensitization. For this reason, the number of channels is allowed to vary as a free parameter in all subsequent analysis. When 23 of the 80 1D-I data sets described above (100 channels per set) are fit to the correct model with the number of channels allowed to vary as a free parameter, the estimate of the number of channels is 100.8 ± 1.1. This suggests that CVF can accurately estimate the number of channels in simulated data generated using a simple model.
Analysis of simulated data suggests that at least 200 points are needed to reliably identify the presence of a fourth desensitized state (see above). For this reason, six GABAR data sets of 200 points are evenly selected from the six raw traces. Reapplication of CVF to these six sets with the number of channels as a ninth free parameter still favors 3D-I,ddp.
When the same six data sets are fit to models with two open states, the most favored model (referred to as COO; shown in Fig. 7) has two open states in a linear arrangement and three desensitized states, one connected to the closed and one connected to each of the two open states. This model is analogous to 3D-I, in that there is one proximal and two distal desensitized states. COO is favored over 3D-I with a median L[LR] of 3.9 and six out of six above 0 (p = 0.016).
|
|
Initially, 200 points are selected evenly along the path of desensitization from the raw simulated COO data. Gaussian noise is added (
= 0.564 channels, average of the baseline noise of the six GABAR traces) and all 10 data sets are fit to model COO. There is a comparable variability in the parameter estimates, particularly regarding
2 (the transitions leading from the second open-state back to the first open state). When the two highest and two lowest estimates are excluded, the average of the central six estimates (42 ± 21 s1) still differs by an unacceptable degree from the simulation parameter (10 s1). To address this issue, three 100-point data sets are selected from the raw COO data with timescale factor (see Methods) values of 1, 3, and 10 and fit to COO. The resulting averages for the central six estimates of
2 are: 37 ± 5 s1, 24 ± 4 s1, and 69 ± 20 s1, respectively. Because the middle value is closest to the simulation parameter, tsf = 3 is chosen for subsequent data point selection. When 200 points are selected from the raw COO data with tsf = 3, the estimate of
2 based on the central six estimates is within an acceptable range of the true value (16 ± 3 s1) and nine of the 11 parameter estimates are within two SEs of the true value (Table 7). The above results on COO data suggest that 200 points selected with tsf = 3 is the optimal method of point selection in the analysis of GABAR data. The final analysis of both GABAR and COO data is done using this point-selection method. Unfortunately, a comparable degree of variation in the parameter estimates for GABAR data remained despite this change (Table 7).
Analysis using optimally selected points
Fig. 6 shows a representative GABAR trace before and after selection of 200 points with tsf = 3, as well as the results of fitting the selected data to the COO model. The COO model is heavily favored over 3D-I, and similar results are observed when fitting simulated COO data but not when fitting simulated 3D-I data (Table 2). When comparing COO to a model containing an additional open state connected to the closed state, two of the six GABAR sets had L[LR] of <0.1. Only two of 10 COO sets had L[LR] of <0.1 when making the same comparison. This differs significantly (p = 0.044) from the prediction that half the L[LR]s are <0.1 when comparing a correct model with a model containing an additional desensitized state. The addition of an unnecessary open state may have a greater effect on L[LR] than the addition of an unnecessary desensitized state. When fitting simulated COO data, the total L[LH] is also greater when an unnecessary open state is added than when an unnecessary desensitized state is added (Fig. 7 B).
For both GABAR and simulated COO data, 2D models are heavily disfavored, whereas at least half of the comparisons involving 4D models have L[LR]s <0.1 (Table 2). This includes a model analogous to the one reported by Haas and Macdonald (1999)
with two nonconducting states separating the two open states (OCDO). Because COO is not a submodel of OCDO, some L[LR]s are <0 (Table 2).
Fig. 7 shows the sum of L[LH]s that result from fitting both GABAR and simulated COO data to 3D-I and 11 models with at least two open states. Similar results are observed for both GABAR and COO data. Among the models with 10 kinetic parameters tested, COO is slightly favored over a model with the closed state between the two open states (based on total L[LH] alone) and heavily favored over a model analogous to 3D-II with two desensitized states connected to the closed state (Table 3). For GABAR data, the ddp pattern of desensitization is slightly favored over the pdd pattern (based on the total L[LH] alone, p = 0.094), and for COO data, the ddp pattern is heavily favored (p < 0.01).
When analyzing data recorded from a biological system, such as a population of ion channels, there is no way to be certain that a particular model is the correct model. The eventual goal is to obtain a model or models adequate for the purpose of understanding biologically relevant functions and pharmacological modulation. It is, therefore, useful to determine if CVF can extract meaningful information when fitting to models known to be simpler than the correct model. L[LR]s are determined for comparing the arrangement of states and the pattern of desensitization for GABAR and COO data assuming only a single open state. In both cases, fitting to 3D-I is heavily favored over 3D-II (Table 3), and the ddp pattern of desensitization is heavily favored over the pdd pattern (Table 5). Despite using an incorrect model, the number of desensitized states on each side of the gating isomerization and their relationship to the phases of desensitization are discernable.
Comparing the SEs to the means, parameter estimates of the six GABAR traces show a great deal of variability despite the use of optimally selected points. This is particularly true of the two parameters leading away from the second open state (
2 and Kf3; see Table 7). The variability does not appear to be a simple reflection of a shallow likelihood surface. Changing either parameter by as little as 10% from the value associated with the MLE (leaving the other nine parameters at their MLE associated values) decreases the L[LH] by at least a factor of 0.05 for all six GABAR traces. This value is well within the resolution of CVF (see Fig. 5 B).
Properties of the COO consensus model
To determine the properties of the COO model, the averages of the central four parameter estimates from fitting the six GABAR data sets are used to generate a final consensus model (Table 7 and Fig. 8). To allow simulation of deactivation, a single agonist binding step connected to the closed state is included. An on-rate of 107 M1 s1 and an off-rate of 500 s1 are used (Scheller and Forman, 2002
), giving an EC50 of
20 µM. A series of smooth simulations are done and the results fit to exponential equations by sum-of-squares minimization.
Desensitization in response to 1 mM GABA proceeds with three distinct phases having time constants of 3.6 ms (62%), 36 ms (14%), 362 ms (23%). These values are comparable to the average values observed when the decaying phase of the six raw GABAR traces are fit by sum-of-squares minimization to a three-exponential equation: 3.8 ms (56%), 64 ms (20%), 730 ms (19%).
When the deactivation phase after a 1-ms application of 1 mM GABA is fit to a double exponential, the time constants are 5 and 140 ms with the fast phase representing 69% of the total deactivation. These values are comparable to those reported by others studying recombinant (Haas and Macdonald, 1999
; Tia et al., 1996
) as well as native receptors (Jones and Westbrook, 1995
; Mozrzymas et al., 1999
; Shen et al., 2000
; Zhu and Vicini, 1997
). The relative magnitude of the slow component increases as the duration of GABA application increases (Fig. 8) as has been described for native receptors (Jones and Westbrook, 1995
).
Recovery from desensitization has been reported to proceed in a biphasic manner with fast time constants ranging between 35 and 150 ms (Bai et al., 1999
; Jones and Westbrook, 1995
; Shen et al., 2000
; Zhu and Vicini, 1997
). In its current form, COO does not predict the rapid phase of recovery from desensitization. This is due to accumulation of receptors in the D1 state during deactivation.
First latency analysis applied to recombinant receptors composed of
1ß1
2 subunits favored a 3D-I,pdd model (Burkat et al., 2001
). The pdd pattern was chosen because high concentrations of GABA often show first latencies of tens of milliseconds. They suggest that these long first latencies are the result of channels briefly visiting the proximal fast desensitized state before opening. In its current form, the COO model does not predict these first latencies.
The COO model predicts mean open times of 0.49 and 2.5 ms, which agree well with 0.3 and 1.92 ms, the two fastest values reported for single-channel analysis of
1ß3
2L receptors recorded at steady state in the presence of 1 mM GABA (Haas and Macdonald, 1999
). The predicted mean closed times of 0.19, 0.87, and 8.5 ms also agree well with the values of 0.2, 1.64, and 8.89 ms reported in the same study.
| DISCUSSION |
|---|
|
|
|---|
2-distribution. In the latter case, Monte Carlo analysis was used to determine the expected distribution of the L[LR]. This study applies similar principals to likelihoods identified by fitting macroscopic nonequilibrium data directly to kinetic models.
When comparing fits to a complex model and a correct submodel, the distribution of 2 x L[LR] observed in this study deviates dramatically from that of a
2-distribution with 2 df. This is most likely due to the restrictions on kinetic parameters. The Akaike information criterion (Akaike, 1974
) is often used for comparing models that differ in complexity when the simple model is not a submodel of the complex one. As with the
2-prediction, this value assumes that parameter estimates are approximately normally distributed around their true values. This requires they be unrestricted in the neighborhood of their true values. When a submodel is correct, the true value of at least one unnecessary parameter of the complex model is either 0 (forward rate constant) or infinite (reverse rate constant). Because kinetic parameters cannot be below 0, this parameter cannot be normally distributed around its true values in either case. This results in a much narrower distribution of L[LR] than would be the case if the unnecessary parameters were unrestricted. In this study, this narrow distribution leads to a useful null hypothesis for ruling out simple models.
In this study, comparisons between simple and complex models are made based on results of fitting to simulated data. For all levels of complexity studied, addition of an unnecessary desensitized state, or an unnecessary transition between states leads to L[LR] distributions with a consistent feature: at least half of the L[LR]s fall below 0.1. Using this feature as a null hypothesis allows the calculation of a p-value for ruling out submodels based on a binomial distribution. An increase in the number of channels and the number of points selected increases the magnitude of L[LR]s when the complex model is correct but has no obvious effect when the submodel is correct (Table 2). From this it appears that when analyzing recorded data, the power of CVF to identify the correct number of desensitized states is primarily limited only by the number of channels in the patch and the number of points selected.
Based on fitting of COO data, the addition of an unnecessary open state may be fundamentally different from the addition of a desensitized state. Only two of 10 L[LR]s are below 0.1 when COO and OCOO are compared. Based on the results of fitting the six GABAR sets to a model with three open states, there is no compelling justification for including a third open state at this time. A third open state, however, is far from being ruled out because only one of the many models with three open states has thus far been tested. In addition, the existence of a third open state is supported by the open time distribution reported from single-channel studies (Macdonald et al., 1989
; Weiss and Magleby, 1989
).
Comparing models with the same number of states and free parameters is made possible by the use of the full covariance matrix. Fitting by wSS, which only takes into account the magnitude of the current at each point in time, fails to discriminate between differences in the arrangement of states or in the pattern of desensitization. In making such comparisons, a binomial distribution is used to generate p-values based on the null hypothesis that L[LR]s are distributed around 0. Although this is adequate for some comparisons, analysis of simulated data reveals several instances in which models known to be incorrect are not ruled out by the L[LR] distribution, despite having a lower total L[LH] than the correct model. Analysis of COO data, for example, excludes the OCO model based only on the total L[LH], not on the L[LR] distribution. This problem may prove surmountable by using a more accurate L[LR] distribution as the null hypothesis. This would require analyzing simulated data generated using OCO to determine the L[LR] distribution associated with comparing OCO and COO models when the former is correct. This distribution can then be used instead of the binomial distribution as a null hypothesis to rule out OCO.
CVF is very sensitive to changes in the number of desensitized states on each side of the gating isomerization. However, it has been impossible to discriminate between linear and branched arrangements that have the same number of desensitized states on each side of the gating isomerization. This comparison is fundamentally different from others made. There is no difference in the total L[LH] and for each data set, L[LR] is close to zero. This finding is based on the analysis of simulated data in which there is a single perturbation. It remains to be seen if different application paradigms with multiple perturbations (such as those testing deactivation or recovery from desensitization) prove more powerful in discriminating between such models (see below).
Parameter estimation
ML theory predicts that when fitting to the correct model, the observed parameter estimates will be approximately normally distributed around the true values. Because this is an asymptotic theory, its predictions are only observed when the number of determinations is adequate. In this study, simulation parameter values including the number of channels are usually within two SEs of the mean for all single open-state data sets examined.
Both COO and GABAR data show a great deal of variability in parameter estimates when fitting to the COO model, despite adequate curvature of the likelihood surface. This suggests that when fitting to models with 11 free parameters, CVF is close to the limit of its accuracy when studying single perturbation data, with 200 channels and 200 points per set. Reasonable parameter estimates are obtained by excluding the outer four of 10 estimates. In applying CVF to recorded data, it may be necessary to analyze >10 replicates to derive accurate parameter estimates. Added power is also expected by analyzing data recorded during multiple perturbations and by simultaneously fitting multiple data sets recorded from the same patch (see below).
Several factors may contribute to the variability in parameter estimates observed when fitting GABAR data to COO. Agonist binding is assumed in this study to be instantaneous, whereas the reported activation rates for 1 mM GABA are less than those for 10 mM GABA (Li and Pearce, 2000
; Lavoie et al., 1997
). In addition to lacking the agonist binding step, the COO model is likely to require additional levels of complexity to fully predict the relevant aspects of GABAR function. Finally, the use of a subunit ratio of (1:1:1.2,
1:ß2:
2) HEK cell transfection may have resulted in a heterogeneous mixture of
1ß2 and
1ß2
2 receptors (Boileau et al., 2002
). By comparing variability within a patch to variability between patches, CVF will be able to address the effect of heterogeneity.
COO consensus model
The main goal of this article is to introduce CVF and report its capabilities and limitations. Application of CVF to GABAR data is primarily done for the purpose of guiding the choice of models to be used for generating simulated data. The consensus model presented here is derived from the results of analyzing only six traces recorded in response to a single perturbation. A more comprehensive model will require a wide range of experimental protocols including deactivation, recovery from desensitization, and application of subsaturating GABA concentrations. Despite this, the COO model reproduces several general features of GABAR function. Because CVF can extract meaningful information from simulated data using models known to be simpler than the correct model, the qualitative features of the COO model are worth discussing briefly.
Although the COO model is derived exclusively from desensitization data, it reproduces the biphasic nature of GABAR deactivation. Failure of the COO model to reproduce the fast phase of recovery from desensitization may be due to the lack of a provision for GABA dissociation from any of the desensitized states. In its present form all activation and deactivation go through the C state, which is connected to the most thermodynamically stable desensitized state. During deactivation, desensitization continues to occur. This can be avoided by allowing GABA to dissociate from the desensitized states, as has been demonstrated for receptors containing the
1ß2
2 subunits (Chang et al., 2002
). Failure to reproduce the prolonged first latencies seen for high GABA may be due to the assumption that all receptors start in the activatable state. If a proportion of receptors are in the D2 state at rest, prolonged first latencies could represent recovery from D2 to O1 after agonist binding. This can be tested by allowing the fraction of receptors in the D2 state at rest to vary as a free parameter.
The kinetic parameters of the COO consensus model predict single-channel properties remarkably similar to those reported for receptors containing the
1ß3
2 subunits in the presence of high GABA at steady state (Haas and Macdonald, 1999
). This same report includes evidence of a third open state and two nonconducting states separating two open states. When tested here, neither feature is statistically justified or definitively ruled out.
Future directions
Although statistical power increases with the number of points selected, this approach is limited by the computational demands of CVF that increase with the square of the number of points selected. The total number of points fit, however, can be increased by simultaneously fitting two traces from the same patch. As long as the two traces are separated by enough time, their covariance matrices can be calculated separately. The likelihood of observing both traces is simply the product of the likelihood of observing the individual traces. Using traces with different paradigms can emphasize different aspects of a complex model. The same is true for more than two traces and the computational demand only goes up linearly with each added trace.
The CVF method will also be expanded to analyze data sets with more than one perturbation. Each perturbation represents a change in the GABA concentration that is represented by the use of a new Q matrix. The new Q matrix is identical to the first except for the elements corresponding to GABA association. The covariance formula (Eq. 4, Appendix I) can be expanded to calculate the covariance between points on different sides of the second perturbation. Equation 4 will be used as is for all pairs of points on the same side of the second perturbation. This can be used to study deactivation as well as any paradigm that involves changing from one GABA concentration to another. A similar process will be used to study three perturbations that will allow analysis of data studying recovery from desensitization. These added paradigms are expected to emphasize different aspects of receptor function. Parameters poorly estimated using one paradigm might be better estimated using another.
Single-channel analysis has long been the most reliable method for obtaining information on molecular mechanism and discriminating between complex kinetic models. CVF reveals that macroscopic traces contain more information than was previously realized and this study begins to explore its capabilities. In addition to complementing single-channel analysis, CVF has several important advantages. The biological importance of ligand-gated ion channels is their role in mediating fast synaptic potentials. Vital to this end is understanding their response to rapid changes in agonist concentration. CVF is ideally suited for studying ion channel kinetics in response to multiple per