Parameterization for In-Silico Modeling of Ion Channel Interactions with Drugs

Since the first Hodgkin and Huxley ion channel model was described in the 1950s, there has been an explosion in mathematical models to describe ion channel function. As experimental data has become richer, models have concomitantly been improved to better represent ion channel kinetic processes, although these improvements have generally resulted in more model complexity and an increase in the number of parameters necessary to populate the models. Models have also been developed to explicitly model drug interactions with ion channels. Recent models of drug-channel interactions account for the discrete kinetics of drug interaction with distinct ion channel state conformations, as it has become clear that such interactions underlie complex emergent kinetics such as use-dependent block. Here, we describe an approach for developing a model for ion channel drug interactions. The method describes the process of extracting rate constants from experimental electrophysiological function data to use as initial conditions for the model parameters. We then describe implementation of a parameter optimization method to refine the model rate constants describing ion channel drug kinetics. The algorithm takes advantage of readily available parallel computing tools to speed up the optimization. Finally, we describe some potential applications of the platform including the potential for gaining fundamental mechanistic insights into ion channel function and applications to in silico drug screening and development.     

Molecular variability elicits a tunable switch with discrete neuromodulatory response phenotypes

Conductance-based compartment modeling requires tuning of many parameters to fit the neuron model to target electrophysiological data. Automated parameter optimization via evolutionary algorithms (EAs) is a common approach to accomplish this task, using error functions to quantify differences between model and target. We present a three-stage EA optimization protocol for tuning ion channel conductances and kinetics in a generic neuron model with minimal manual intervention. We use the technique of Latin hypercube sampling in a new way, to choose weights for error functions automatically so that each function influences the parameter search to a similar degree. This protocol requires no specialized physiological data collection and is applicable to commonly-collected current clamp data and either single- or multi-objective optimization. We applied the protocol to two representative pyramidal neurons from layer 3 of the prefrontal cortex of rhesus monkeys, in which action potential firing rates are significantly higher in aged compared to young animals. Using an idealized dendritic topology and models with either 4 or 8 ion channels (10 or 23 free parameters respectively), we produced populations of parameter combinations fitting the target datasets in less than 80 hours of optimization each. Passive parameter differences between young and aged models were consistent with our prior results using simpler models and hand tuning. We analyzed parameter values among fits to a single neuron to facilitate refinement of the underlying model, and across fits to multiple neurons to show how our protocol will lead to predictions of parameter differences with aging in these neurons.     

Gating kinetics of potassium channel in rat dorsal root ganglion neurons analyzed with fractal model

The kinetics of ion channels have been widely modeled as a Markov process. In these models it is assumed that the channel protein has a small number of discrete conformational states and kinetic rate constants connecting these states are constant. To study the gating kinetics of voltage-dependent K(+) channel in rat dorsal root ganglion neurons, K(+) channel current were recorded using cell-attached patch-clamp technique. The K(+) channel characteristic of kinetics were found to be statistically self-similar at different time scales as predicted by the fractal model. The fractal dimension D for the closed times and for the open times depend on the pipette potential. For the open and closed times of kinetic setpoint, it was found dependent on the applied pipette potential, which indicated that the ion channel gating kinetics had nonlinear kinetic properties. Thus, the open and closed durations, which had the voltage dependence of the gating of this ion channel, were well described by the fractal model.     

A numerical approach to ion channel modelling using whole-cell voltage-clamp recordings and a genetic algorithm

The activity of trans-membrane proteins such as ion channels is the essence of neuronal transmission. The currently most accurate method for determining ion channel kinetic mechanisms is single-channel recording and analysis. Yet, the limitations and complexities in interpreting single-channel recordings discourage many physiologists from using them. Here we show that a genetic search algorithm in combination with a gradient descent algorithm can be used to fit whole-cell voltage-clamp data to kinetic models with a high degree of accuracy. Previously, ion channel stimulation traces were analyzed one at a time, the results of these analyses being combined to produce a picture of channel kinetics. Here the entire set of traces from all stimulation protocols are analysed simultaneously. The algorithm was initially tested on simulated current traces produced by several Hodgkin-Huxley–like and Markov chain models of voltage-gated potassium and sodium channels. Currents were also produced by simulating levels of noise expected from actual patch recordings. Finally, the algorithm was used for finding the kinetic parameters of several voltage-gated sodium and potassium channels models by matching its results to data recorded from layer 5 pyramidal neurons of the rat cortex in the nucleated outside-out patch configuration. The minimization scheme gives electrophysiologists a tool for reproducing and simulating voltage-gated ion channel kinetics at the cellular level.     

Fractal analysis of a voltage-dependent potassium channel from cultured mouse hippocampal neurons.

The kinetics of ion channels have been widely modeled as a Markov process. In these models it is assumed that the channel protein has a small number of discrete conformational states and the kinetic rate constants connecting these states are constant. In the alternative fractal model the spontaneous fluctuations of the channel protein at many different time scales are represented by a kinetic rate constant k = At1-D, where A is the kinetic setpoint and D the fractal dimension. Single-channel currents were recorded at 146 mM external K+ from an inwardly rectifying, 120 pS, K+ selective, voltage-sensitive channel in cultured mouse hippocampal neurons. The kinetics of these channels were found to be statistically self-similar at different time scales as predicted by the fractal model. The fractal dimensions were approximately 2 for the closed times and approximately 1 for the open times and did not depend on voltage. For both the open and closed times the logarithm of the kinetic setpoint was found to be proportional to the applied voltage, which indicates that the gating of this channel involves the net inward movement of approximately one negative charge when this channel opens. Thus, the open and closed times and the voltage dependence of the gating of this channel are well described by the fractal model.     

A model for conductance changes in the squid giant axon. II. Channel kinetics

The kinetics of the sodium and potassium channels in voltage clamped squid giant axon following a relaxation of the membrane subunits are examined and compared with the Hodgkin-Huxley equations. Mechanisms are suggested for the turn-off of the sodium conductance and a set of kinetic states are proposed for the potassium channel which are consistent with the experimental observations. Determination of the rate constants for relaxation of the surface subunits which triggers the subsequent changes within the independent channels provide information on the equilibrium constant and free energy for this process. The free energy is observed to approach zero as the depolarizing voltage of the clamp approaches $\text{E}{}_{\mathrm{<mtext>Na</mtext>}}$, the voltage for zero sodium current in voltage clamped axons. Analysis of the final rate constants in the kinetic sequence for potassium indicates a symmetry of the channel when it is in its steady-state configuration during clamp in the absence of external gradients.     

Relationship between membrane excitability and single channel open-close kinetics.

We have developed a novel technique for simulating the influence of the effects of single channel kinetics on the voltage changes associated with membrane excitability. The technique uses probability distribution functions for the durations of channel open- and closed-state lifetimes, which can be calculated for any model of the ion conductance process. To illustrate the technique, we have used the Hodgkin and Huxley model of nerve membrane ion conductances to simulate channel kinetics during predetermined voltage changes, such as a voltage jump and an action potential. We have also simulated the influence of channels on voltage changes in a free running, non-voltage-clamped patch of membrane 1 micron2 or less in area. The latter results provide a direct illustration of the relationship between fluctuations of membrane excitability and fluctuations in channel open- and closed-state lifetimes.     

Exploring voltage-dependent ion channels in silico by hysteretic conductance

Kinetic models of voltage-dependent ion channels are normally inferred from time records of macroscopic current relaxation or microscopic single channel data. A complementary explorative approach is outlined. Hysteretic conductance refers to conductance delays in response to voltage changes, delays at either macroscopic or microscopic levels of observation. It enables complementary assessments of model assumptions and gating schemes of voltage-dependent channels, e.g. independent versus cooperative gating, and multiple gating modes. Under the Hodgkin-Huxley condition of independent gating, and under ideal measurement conditions, hysteretic conductance makes it also possible to estimate voltage-dependent rate functions. The argument is mainly theoretical, based on experimental observations, and illustrated by simulations of Markov kinetic models.     

Dynamical estimation of neuron and network properties II: path integral Monte Carlo methods

Hodgkin–Huxley (HH) models of neuronal membrane dynamics consist of a set of nonlinear differential equations that describe the time-varying conductance of various ion channels. Using observations of voltage alone we show how to estimate the unknown parameters and unobserved state variables of an HH model in the expected circumstance that the measurements are noisy, the model has errors, and the state of the neuron is not known when observations commence. The joint probability distribution of the observed membrane voltage and the unobserved state variables and parameters of these models is a path integral through the model state space. The solution to this integral allows estimation of the parameters and thus a characterization of many biological properties of interest, including channel complement and density, that give rise to a neuron’s electrophysiological behavior. This paper describes a method for directly evaluating the path integral using a Monte Carlo numerical approach. This provides estimates not only of the expected values of model parameters but also of their posterior uncertainty. Using test data simulated from neuronal models comprising several common channels, we show that short (<50 ms) intracellular recordings from neurons stimulated with a complex time-varying current yield accurate and precise estimates of the model parameters as well as accurate predictions of the future behavior of the neuron. We also show that this method is robust to errors in model specification, supporting model development for biological preparations in which the channel expression and other biophysical properties of the neurons are not fully known.     

Electrical properties of ionic channels formed by Helix pomatia hemocyanin in planar lipid bilayers

Helix pomatia hemocyanin forms ion-conducting channels in planar lipid bilayer membranes when added at mg/ml concentration. These channels have several original features. They fluctuate between one conducting and some poorly conducting states and fluctuations can be grouped in bursts. Different channels can have widely different conductance amplitudes. Both channel conductance and burst lifetime are dependent on the applied voltage. Fluctuations within a burst show a complex kinetic behaviour which has been explained developing a multistate model. The model calls for one single open state and six different closed states. Transitions are allowed only between one of the closed states and the open one and obey first order kinetics. This model is able to fit all our experimental curves obtained in single channel experiments.     

Estimation of ion channel kinetics from fluctuations of macroscopic currents

For single channel recordings, the maximum likelihood estimation (MLE) of kinetic rates and conductance is well established. A direct extrapolation of this method to macroscopic currents is computationally prohibitive: it scales as a power of the number of channels. An approximated MLE that ignored the local time correlation of the data has been shown to provide estimates of the kinetic parameters. In this article, an improved approximated MLE that takes into account the local time correlation is proposed. This method estimates the channel kinetics using both the time course and the random fluctuations of the macroscopic current generated by a homogeneous population of ion channels under white noise. It allows arbitrary kinetic models and stimulation protocols. The application of the proposed algorithm to simulated data from a simple three-state model on nonstationary conditions showed reliable estimates of all the kinetic constants, the conductance and the number of channels, and reliable values for the standard error of those estimates. Compared to the previous approximated MLE, it reduces by a factor of 10 the amount of data needed to secure a given accuracy and it can even determine the kinetic rates in macroscopic stationary conditions.     

Properties and modulation of alpha human atrial natriuretic peptide (alpha-hANP)-formed ion channels

Using the lipid bilayer technique we have optimized recording conditions and confirmed that alpha human atrial natriuretic peptide [alpha-hANP(1-28)] forms single ion channels. The single channel currents recorded in 250/50 mM KCl cis/trans chambers show that the ANP-formed channels were heterogeneous, and differed in their conductance, kinetic, and pharmacological properties. The ANP-formed single channels were grouped as: (i) H202- and Ba2+-sensitive channel with fast kinetics; the nonlinear current-voltage (I-V) relationship of this channel had a reversal potential (Erev) of -28.2 mV, which is close to the equilibrium potential for K+ (EK = -35 mV) and a maximal slope conductance (gmax) of 68 pS at positive potentials. Sequential ionic substitution (KCl, K gluconate and choline Cl) of the cis solution suggests that the current was carried by cations. The fast channel had three modes (spike mode, burst mode, and open mode) that differed in their kinetics but not in their conductance properties. (ii) A large conductance channel possessing several subconductance levels that showed time-dependent inactivation at positive and negative membrane potentials (Vm). The inactivation ratio of the current at the end of the voltage step (Iss) to the initial current (Ii) activated immediately after the voltage step, (Iss/Ii), was voltage dependent and described by a bell-shaped curve. The maximal current-voltage (I-V) relationship of this channel, which had an Erev of +17.2 mV, was nonlinear and the value of gmax was 273 pS at negative voltages. (iii) A transiently-activated channel: the nonlinear I-V relationship of this channel had an Erev of -29.8 mV and the value of gmax was 160 pS at positive voltages. We propose that the voltage-dependence of the ionic currents and the kinetic parameters of these channel types indicate that if they were formed in vivo and activated by cytosolic factors they could change the membrane potential and the electrolyte homeostasis of the cell.     

Voltage-dependence of ion permeation in cyclic GMP-gated ion channels is optimized for cell function in rod and cone photoreceptors

The kinetics of the photocurrent in both rod and cone retinal photoreceptors are independent of membrane voltage over the physiological range (-30 to -65 mV). This is surprising since the photocurrent time course is regulated by the influx of Ca(2+) through cGMP-gated ion channels (CNG) and the force driving this flux changes with membrane voltage. To understand this paradigm, we measured Pf, the fraction of the cyclic nucleotide-gated current specifically carried by Ca(2+) in intact, isolated photoreceptors. To measure Pf we activated CNG channels by suddenly increasing free 8-Br-cGMP in the cytoplasm of rods or cones loaded with a caged ester of the cyclic nucleotide. Simultaneous with the uncaging flash, we measured the cyclic nucleotide-dependent changes in membrane current and fluorescence of the Ca(2+) binding dye, Fura-2, also loaded into the cells. We determined Pf under physiological solutions at various holding membrane voltages between -65 and -25 mV. Pf is larger in cones than in rods, but in both photoreceptor types its value is independent of membrane voltage over the range tested. This biophysical feature of the CNG channels offers a functional advantage since it insures that the kinetics of the phototransduction current are controlled by light, and not by membrane voltage. To explain our observation, we developed a rate theory model of ion permeation through CNG channels that assumes the existence of two ion binding sites within the permeation pore. To assign values to the kinetic rates in the model, we measured experimental I-V curves in membrane patches of rods and cones over the voltage range -90 to 90 mV in the presence of simple biionic solutions at different concentrations. We optimized the fit between simulated and experimental data. Model simulations describe well experimental photocurrents measured under physiological solutions in intact cones and are consistent with the voltage-independence of Pf, a feature that is optimized for the function of the channel in photoreceptors.     

Applications of nonequilibrium response spectroscopy to the study of channel gating. Experimental design and optimization

A novel experimental technique known as non-equilibrium response spectroscopy (NRS) based on ion channel responses to rapidly fluctuating voltage waveforms was recently described (Millonas & Hanck, 1998a). It was demonstrated that such responses can be affected by subtle details of the kinetics that are otherwise invisible when conventional stepped pulses are applied. As a consequence, the kinetics can be probed in a much more sensitive way by supplementing conventional techniques with measurements of the responses to more complex voltage waveforms. In this paper we provide an analysis of the problem of the design and optimization of such waveforms. We introduce some methods for determination of the parametric uncertainty of a class of kinetic models for a particular data set. The parametric uncertainty allows for a characterization of the amount of kinetic information acquired through a set of experiments which can in turn be used to design new experiments that increase this information. We revisit the application of dichotomous noise (Millonas & Hanck, 1998a, b), and further consider applications of a more general class of continuous wavelet -based waveforms. A controlled illustration of these methods is provided by making use of a simplified "toy" model for the potassium channel kinetics.     

Efficient fitting of conductance-based model neurons from somatic current clamp

Estimating biologically realistic model neurons from electrophysiological data is a key issue in neuroscience that is central to understanding neuronal function and network behavior. However, directly fitting detailed Hodgkin?CHuxley type model neurons to somatic membrane potential data is a notoriously difficult optimization problem that can require hours/days of supercomputing time. Here we extend an efficient technique that indirectly matches neuronal currents derived from somatic membrane potential data to two-compartment model neurons with passive dendrites. In consequence, this approach can fit semi-realistic detailed model neurons in a few minutes. For validation, fits are obtained to model-derived data for various thalamo-cortical neuron types, including fast/regular spiking and bursting neurons. A key aspect of the validation is sensitivity testing to perturbations arising in experimental data, including sampling rates, inadequately estimated membrane dynamics/channel kinetics and intrinsic noise. We find that maximal conductance estimates and the resulting membrane potential fits diverge smoothly and monotonically from near-perfect matches when unperturbed. Curiously, some perturbations have little effect on the error because they are compensated by the fitted maximal conductances. Therefore, the extended current-based technique applies well under moderately inaccurate model assumptions, as required for application to experimental data. Furthermore, the accompanying perturbation analysis gives insights into neuronal homeostasis, whereby tuning intrinsic neuronal properties can compensate changes from development or neurodegeneration.     

Optimal estimation of ion-channel kinetics from macroscopic currents

Markov modeling provides an effective approach for modeling ion channel kinetics. There are several search algorithms for global fitting of macroscopic or single-channel currents across different experimental conditions. Here we present a particle swarm optimization(PSO)-based approach which, when used in combination with golden section search (GSS), can fit macroscopic voltage responses with a high degree of accuracy (errors within 1%) and reasonable amount of calculation time (less than 10 hours for 20 free parameters) on a desktop computer. We also describe a method for initial value estimation of the model parameters, which appears to favor identification of global optimum and can further reduce the computational cost. The PSO-GSS algorithm is applicable for kinetic models of arbitrary topology and size and compatible with common stimulation protocols, which provides a convenient approach for establishing kinetic models at the macroscopic level.     

Corrections for space-clamp errors in measured parameters of voltage-dependent conductances in a cylindrical neurite

The ability to correct parameters of voltage-gated conductances measured under poor spatial control by point voltage clamp could rescue much flawed experimental data. We explore a strategy for correcting errors in experiments that employs a full-trace approach to parameter determination. Simulated soma voltage-clamp runs are made on a model neuron with a single voltage-gated, Hodgkin-Huxley channel type distributed uniformly along an elongate process. Estimates for both kinetic and I(V) parameters are obtained by fitting a form of the Hodgkin-Huxley equations to the complete time course of leak-subtracted current curves. The fitted parameters are used to determine how much correction in each parameter is needed to regenerate the set actually belonging to the channel. Corrections are generated for a range of neurite lengths, conductance densities, and channel characteristics.     

Modulation of gramicidin A open channel lifetime by ion occupancy.

The hypothesis that the gramicidin A channel stability depends on the level of ion occupancy of the channel was used to derive a mathematical model relating channel lifetime to channel occupancy. Eyring barrier permeation models were examined for their ability to fit the zero-voltage conductance, current-voltage, as well as lifetime data. The simplest permeation model required to explain the major features of the experimental data consists of three barriers and four sites (3B4S) with a maximum of two ions occupying the channel. The average lifetime of the channel was calculated from the barrier model by assuming the closing rate constant to be proportional to the probability of the internal channel sites being empty. The link between permeation and lifetime has as its single parameter the experimentally determined averaged lifetime of gramicidin A channels in the limit of infinitely dilute solutions and has therefore no adjustable parameters. This simple assumption that one or more ions inside the channel completely stabilize the dimer conformation is successful in explaining the experimental data considering the fact that this model for stabilization is independent of ion species and configurational occupancy. The model is used to examine, by comparison with experimental data, the asymmetrical voltage dependence of the lifetime in asymmetrical solutions, the effects of blockers, and the effects of elevated osmotic pressure.     

Application of the ensemble nonequilibrium response spectroscopy to shaker potassium channel gating

Standard electrophysiology techniques study relaxation transients in voltage-gated ion channels generated by discrete voltage steps. The nonequilibrium response spectroscopy involves analyzing responses to fluctuating potentials. We apply the ensemble NRS method to gating kinetics of Shaker potassium ion channels. We evaluate various proposed Markov models of channel gating from the nonequilibrium response viewpoint. These new NRS protocols can be used to test otherwise indistinguishable models or improve estimates for parameters of channel kinetics models.     

A multi-ion permeation mechanism in neuronal background chloride channels

Unitary current/voltage relationships of background Cl channels of rat hippocampal neurons were determined for varied gradients and absolute concentrations of NaCl. The channels revealed permeabilities for both Cl and Na ions. A hyperlinear increase of unitary conductance, observed for a symmetrical increase of salt concentration from 300 and 600 mM, indicated a multi-ion permeation mechanism. A variety of kinetic models of permeation were tested against the experimental current/voltage relationships. Models involving a pore occupied by mixed complexes of up to five ions were necessary to reproduce all measurements. A minimal model included four equilibrium states and four rate-limiting transitions, such that the empty pore accepts first an anion and then can acquire one or two cation/anion pairs. Three transport cycles are formed: a slow anion cycle (between the empty and single-anion states), a slow cation cycle (between the one- and three-ion states), and a fast anion cycle (between the three- and five-ion states). Thus, permeant anions are required for cation permeation, and several bound anions and cations promote a high rate of anion permeation. The optimized free- energy and electrical charge parameters yielded a self-consistent molecular interpretation, which can account for the particular order in which the pore accepts ions from the solutions. Although the model describes the mixed anion/cation permeability of the channel observed at elevated concentrations, it predicts a high selectivity for Cl anion at physiological ionic conditions.