A general approach to modeling conduction and concentration dynamics in excitable cells of concentric cylindrical geometry

Membrane electrical properties of mouse hippocampal CA1 pyramidal neurons during strong inputs

In this work, we highlight an electrophysiological feature often observed in recordings from mouse CA1 pyramidal cells that has so far been ignored by experimentalists and modelers. It consists of a large and dynamic increase in the depolarization baseline (i.e., the minimum value of the membrane potential between successive action potentials during a sustained input) in response to strong somatic current injections. Such an increase can directly affect neurotransmitter release properties and, more generally, the efficacy of synaptic transmission. However, it cannot be explained by any currently available conductance-based computational model. Here we present a model addressing this issue, demonstrating that experimental recordings can be reproduced by assuming that an input current modifies, in a time-dependent manner, the electrical and permeability properties of the neuron membrane by shifting the ionic reversal potentials and channel kinetics. For this reason, we propose that any detailed model of ion channel kinetics for neurons exhibiting this characteristic should be adapted to correctly represent the response and the synaptic integration process during strong and sustained inputs.

3D Ultrastructure of the Cochlear Outer Hair Cell Lateral Wall Revealed By Electron Tomography

Outer Hair Cells (OHCs) in the mammalian cochlea display a unique type of voltage-induced mechanical movement termed electromotility, which amplifies auditory signals and contributes to the sensitivity and frequency selectivity of mammalian hearing. Electromotility occurs in the OHC lateral wall, but it is not fully understood how the supramolecular architecture of the lateral wall enables this unique form of cellular motility. Employing electron tomography of high-pressure frozen and freeze-substituted OHCs, we visualized the 3D structure and organization of the membrane and cytoskeletal components of the OHC lateral wall. The subsurface cisterna (SSC) is a highly prominent feature, and we report that the SSC membranes and lumen possess hexagonally ordered arrays of particles. We also find the SSC is tightly connected to adjacent actin filaments by short filamentous protein connections. Pillar proteins that join the plasma membrane to the cytoskeleton appear as variable structures considerably thinner than actin filaments and significantly more flexible than actin-SSC links. The structurally rich organization and rigidity of the SSC coupled with apparently weaker mechanical connections between the plasma membrane (PM) and cytoskeleton reveal that the membrane-cytoskeletal architecture of the OHC lateral wall is more complex than previously appreciated. These observations are important for our understanding of OHC mechanics and need to be considered in computational models of OHC electromotility that incorporate subcellular features.

Linking demyelination to compound action potential dispersion with a spike-diffuse-spike approach

To establish and exploit novel biomarkers of demyelinating diseases requires a mechanistic understanding of axonal propagation. Here, we present a novel computational framework called the stochastic spike-diffuse-spike (SSDS) model for assessing the effects of demyelination on axonal transmission. It models transmission through nodal and internodal compartments with two types of operations: a stochastic integrate-and-fire operation captures nodal excitability and a linear filtering operation describes internodal propagation. The effects of demyelinated segments on the probability of transmission, transmission delay and spike time jitter are explored. We argue that demyelination-induced impedance mismatch prevents propagation mostly when the action potential leaves a demyelinated region, not when it enters a demyelinated region. In addition, we model sodium channel remodeling as a homeostatic control of nodal excitability. We find that the effects of mild demyelination on transmission probability and delay can be largely counterbalanced by an increase in excitability at the nodes surrounding the demyelination. The spike timing jitter, however, reflects the level of demyelination whether excitability is fixed or is allowed to change in compensation. This jitter can accumulate over long axons and leads to a broadening of the compound action potential, linking microscopic defects to a mesoscopic observable. Our findings articulate why action potential jitter and compound action potential dispersion can serve as potential markers of weak and sporadic demyelination.

Minimizing the caliber of myelinated axons by means of nodal constrictions

In myelinated axons, most of the voltage-gated ion channels are concentrated at the nodes of Ranvier, which are short gaps in the myelin sheath. This arrangement leads to saltatory conduction and a larger conduction velocity than in nonmyelinated axons. Intriguingly, axons in the peripheral nervous system that exceed about 2 μm in diameter exhibit a characteristic narrowing of the axon at nodes that results in a local reduction of the axonal cross-sectional area. The extent of constriction increases with increasing internodal axonal caliber, reaching a threefold reduction in diameter for the largest axons. In this paper, we use computational modeling to investigate the effect of nodal constrictions on axonal conduction velocity. For a fixed number of ion channels, we find that there is an optimal extent of nodal constriction which minimizes the internodal axon caliber that is required to achieve a given target conduction velocity, and we show that this is sensitive to the precise geometry of the axon and myelin sheath in the flanking paranodal regions. Thus axonal constrictions at nodes of Ranvier appear to be a biological adaptation to minimize axonal volume, thereby maximizing the spatial and metabolic efficiency of these processes, which can be a significant evolutionary constraint. We show that the optimal nodal morphologies are relatively insensitive to changes in the number of nodal sodium channels.

Electrodiffusive model for astrocytic and neuronal ion concentration dynamics

The cable equation is a proper framework for modeling electrical neural signalling that takes place at a timescale at which the ionic concentrations vary little. However, in neural tissue there are also key dynamic processes that occur at longer timescales. For example, endured periods of intense neural signaling may cause the local extracellular K⁺-concentration to increase by several millimolars. The clearance of this excess K⁺ depends partly on diffusion in the extracellular space, partly on local uptake by astrocytes, and partly on intracellular transport (spatial buffering) within astrocytes. These processes, that take place at the time scale of seconds, demand a mathematical description able to account for the spatiotemporal variations in ion concentrations as well as the subsequent effects of these variations on the membrane potential. Here, we present a general electrodiffusive formalism for modeling of ion concentration dynamics in a one-dimensional geometry, including both the intra- and extracellular domains. Based on the Nernst-Planck equations, this formalism ensures that the membrane potential and ion concentrations are in consistency, it ensures global particle/charge conservation and it accounts for diffusion and concentration dependent variations in resistivity. We apply the formalism to a model of astrocytes exchanging ions with the extracellular space. The simulations show that K⁺-removal from high-concentration regions is driven by a local depolarization of the astrocyte membrane, which concertedly (i) increases the local astrocytic uptake of K⁺, (ii) suppresses extracellular transport of K⁺, (iii) increases axial transport of K⁺ within astrocytes, and (iv) facilitates astrocytic relase of K⁺ in regions where the extracellular concentration is low. Together, these mechanisms seem to provide a robust regulatory scheme for shielding the extracellular space from excess K⁺.

When neurons generate electrical signals they release potassium ions (K⁺) into the extracellular space. During periods of intense neural activity, the local extracellular K⁺ may increase drastically. If it becomes too high, it can lead to neural dysfunction. Astrocytes (a kind of glial cells) are involved in preventing this from happening. Astrocytes can take up excess K⁺, transport it intracellularly, and release it in regions where the concentration is lower. This process is called spatial buffering, and a full mechanistic understanding of it is currently lacking. The aim of this work is twofold: First, we develop a formalism for modeling ion concentration dynamics in the intra- and extracellular space. The formalism is general, and could be used to simulate many cellular processes. It accounts for ion transports due to diffusion (along concentration gradients) as well as electrical migration (along voltage gradients). It extends previous, related formalisms, which have focused only on intracellular dynamics. Secondly, we apply the formalism to model how astrocytes exchange ions with the extracellular space. We conclude that the membrane mechanisms possessed by astrocytes seem optimal for shielding the extracellular space from excess K⁺, and provide a full mechanistic description of the spatial (K⁺) buffering process.

Subtle paranodal injury slows impulse conduction in a mathematical model of myelinated axons

This study explores in detail the functional consequences of subtle retraction and detachment of myelin around the nodes of Ranvier following mild-to-moderate crush or stretch mediated injury. An equivalent electrical circuit model for a series of equally spaced nodes of Ranvier was created incorporating extracellular and axonal resistances, paranodal resistances, nodal capacitances, time varying sodium and potassium currents, and realistic resting and threshold membrane potentials in a myelinated axon segment of 21 successive nodes. Differential equations describing membrane potentials at each nodal region were solved numerically. Subtle injury was simulated by increasing the width of exposed nodal membrane in nodes 8 through 20 of the model. Such injury diminishes action potential amplitude and slows conduction velocity from 19.1 m/sec in the normal region to 7.8 m/sec in the crushed region. Detachment of paranodal myelin, exposing juxtaparanodal potassium channels, decreases conduction velocity further to 6.6 m/sec, an effect that is partially reversible with potassium ion channel blockade. Conduction velocity decreases as node width increases or as paranodal resistance falls. The calculated changes in conduction velocity with subtle paranodal injury agree with experimental observations. Nodes of Ranvier are highly effective but somewhat fragile devices for increasing nerve conduction velocity and decreasing reaction time in vertebrate animals. Their fundamental design limitation is that even small mechanical retractions of myelin from very narrow nodes or slight loosening of paranodal myelin, which are difficult to notice at the light microscopic level of observation, can cause large changes in myelinated nerve conduction velocity.

Modeling for cardiac excitation propagation based on the Nernst-Planck equation and homogenization

The bidomain model is a commonly used mathematical model of the electrical properties of the cardiac muscle that takes into account the anisotropy of both the intracellular and extracellular spaces. However, the equations contain self-contradiction such that the update of ion concentrations does not consider intracellular or extracellular ion movements due to the gradient of electric potential and the membrane charge as capacitive currents in spite of the fact that those currents are taken into account in forming Kirchhoff's first law. To overcome this problem, we start with the Nernst-Planck equation, the ionic conservation law, and the electroneutrality condition at the cellular level, and by introducing a homogenization method and assuming uniformity of variables at the microscopic scale, we derive rational bidomain equations at the macroscopic level.

Regulation of ion gradients across myocardial ischemic border zones: a biophysical modelling analysis

The myocardial ischemic border zone is associated with the initiation and sustenance of arrhythmias. The profile of ionic concentrations across the border zone play a significant role in determining cellular electrophysiology and conductivity, yet their spatial-temporal evolution and regulation are not well understood. To investigate the changes in ion concentrations that regulate cellular electrophysiology, a mathematical model of ion movement in the intra and extracellular space in the presence of ionic, potential and material property heterogeneities was developed. The model simulates the spatial and temporal evolution of concentrations of potassium, sodium, chloride, calcium, hydrogen and bicarbonate ions and carbon dioxide across an ischemic border zone. Ischemia was simulated by sodium-potassium pump inhibition, potassium channel activation and respiratory and metabolic acidosis. The model predicted significant disparities in the width of the border zone for each ionic species, with intracellular sodium and extracellular potassium having discordant gradients, facilitating multiple gradients in cellular properties across the border zone. Extracellular potassium was found to have the largest border zone and this was attributed to the voltage dependence of the potassium channels. The model also predicted the efflux of [Formula: see text] from the ischemic region due to electrogenic drift and diffusion within the intra and extracellular space, respectively, which contributed to [Formula: see text] depletion in the ischemic region.

Automaticity in acute ischemia: bifurcation analysis of a human ventricular model

Acute ischemia (restriction in blood supply to part of the heart as a result of myocardial infarction) induces major changes in the electrophysiological properties of the ventricular tissue. Extracellular potassium concentration ([K(o)(+)]) increases in the ischemic zone, leading to an elevation of the resting membrane potential that creates an "injury current" (I(S)) between the infarcted and the healthy zone. In addition, the lack of oxygen impairs the metabolic activity of the myocytes and decreases ATP production, thereby affecting ATP-sensitive potassium channels (I(Katp)). Frequent complications of myocardial infarction are tachycardia, fibrillation, and sudden cardiac death, but the mechanisms underlying their initiation are still debated. One hypothesis is that these arrhythmias may be triggered by abnormal automaticity. We investigated the effect of ischemia on myocyte automaticity by performing a comprehensive bifurcation analysis (fixed points, cycles, and their stability) of a human ventricular myocyte model [K. H. W. J. ten Tusscher and A. V. Panfilov, Am. J. Physiol. Heart Circ. Physiol. 291, H1088 (2006)] as a function of three ischemia-relevant parameters [K(o)(+)], I(S), and I(Katp). In this single-cell model, we found that automatic activity was possible only in the presence of an injury current. Changes in [K(o)(+)] and I(Katp) significantly altered the bifurcation structure of I(S), including the occurrence of early-after depolarization. The results provide a sound basis for studying higher-dimensional tissue structures representing an ischemic heart.

Persistent reflection underlies ectopic activity in multiple sclerosis: a numerical study

Ectopic activity in multiple sclerosis (MS) patients has been traditionally attributed to hyperexcitability of the demyelinated axon segments. Here, we propose that the same outcome may be the result of persistent reflection--the continuous reactivation of the axonal nodes that limit a demyelinated internodal segment. Using computer simulations, we studied the patterns of impulse propagation for a wide range of electrophysiological conditions. In uniformly myelinated fibers, increasing the temperature enabled successful propagation with no blocks in more severe conditions of demyelination. Secondary activations that were originated at the paranodes were formed for temperatures lower than T = 305 K, and at the condition of high sodium channel excitability. Non-sustained and persistent reflections appeared in the case of focally demyelinated fibers, and only within a narrow range of parameters of high temperature and membrane excitability. Persistent reflection reached steady state in ionic currents within 4 ms, and was characterized with a very high activation frequency of 1.504 (+/- 0.039 kHz. We conclude that persistent reflection is a possible mechanism for ectopic activity in MS patients, being more prominent in higher temperatures and severe axonal demyelination. Eliminating these symptoms may be addressed by cooling the body or by applying pharmacological agents to alter excitability properties.

Simulations of reduced conduction reserve in the diabetic rat heart: response to uncoupling and reduced excitability

Experimental results have shown that action potential (AP) conduction in ventricular tissue from streptozotocin-diabetic (STZ) rats is compromised. This was manifest as increased sensitivity of conduction velocity (CV) to the gap junction uncoupler heptanol, as well as increased sensitivity of CV to reduced cellular excitability due to elevated extracellular K(+) concentration, in the STZ hearts. This "reduced conduction reserve" has been suggested to be due to lateralization of connexin43 (Cx43) proteins, rendering them nonfunctional, resulting in compromised intercellular electrical coupling. In this study, we have used computer simulations of one-dimensional AP conduction in a model of rat ventricular myocytes to verify this interpretation. Our results show that compromised intercellular coupling indeed reduces conduction reserve and predict a response to gap junction uncoupling with heptanol that is consistent with experiments. However, our simulations also show that compromised intercellular coupling is insufficient to explain the increased sensitivity to reduced cellular excitability. A thorough investigation of possible underlying mechanisms, suggests that subtle alterations in the voltage-dependence of steady-state gating for the Na(+) current (I (Na)), combined with compromised intercellular coupling, is a likely mechanism for these observations.

Transmission-line model for myelinated nerve fiber

Herein, the well-known cable equation for non-myelinated axon model is extended analytically for myelinated axon formulation. The classical cable equation is thereby modified into a linear second order ordinary differential equation with periodic coefficient, known as Hill's equation. Hill's equation exhibits periodic solutions, known as Floquet's modes. The Floquet's modes are recognized as the nerve fiber activation modes, which are conventionally associated with the nonlinear Hodgkin-Huxley formulation. They can also be incorporated in our linear model.

Generalized Cable Equation Model for Myelinated Nerve Fiber

Herein, the well-known cable equation for nonmyelinated axon model is extended analytically for myelinated axon formulation. The myelinated membrane conductivity is represented via the Fourier series expansion. The classical cable equation is thereby modified into a linear second order ordinary differential equation with periodic coefficients, known as Hill's equation. The general internal source response, expressed via repeated convolutions, uniformly converges provided that the entire periodic membrane is passive. The solution can be interpreted as an extended source response in an equivalent nonmyelinated axon (i.e., the response is governed by the classical cable equation). The extended source consists of the original source and a novel activation function, replacing the periodic membrane in the myelinated axon model. Hill's equation is explicitly integrated for the specific choice of piecewise constant membrane conductivity profile, thereby resulting in an explicit closed form expression for the transmembrane potential in terms of trigonometric functions. The Floquet's modes are recognized as the nerve fiber activation modes, which are conventionally associated with the nonlinear Hodgkin-Huxley formulation. They can also be incorporated in our linear model, provided that the periodic membrane point-wise passivity constraint is properly modified. Indeed, the modified condition, enforcing the periodic membrane passivity constraint on the average conductivity only leads, for the first time, to the inclusion of the nerve fiber activation modes in our novel model. The validity of the generalized transmission-line and cable equation models for a myelinated nerve fiber, is verified herein through a rigorous Green's function formulation and numerical simulations for transmembrane potential induced in three-dimensional myelinated cylindrical cell. It is shown that the dominant pole contribution of the exact modal expansion is the transmembrane potential solution of our generalized model.

Internodal sodium channels ensure active processes under myelin manifesting in depolarizing afterpotentials

The current opinion about processes in myelinated axon is that action potential saltatorially propagates between nodes of Ranvier and passively charges internodal axolemma thus causing depolarizing afterpotentials (DAP). Demyelination blocks the conduction that gives additional argument in favor of hypothesis that internode is not able to be activated by the existing internodal sodium channels. The results of our modeling study shows that, when periaxonal space is sufficiently narrow, saltatorial action potential is able to activate internodes. Low density of internodal sodium channels is sufficient to generate active internodal waves that slowly propagate from nodes towards corresponding midinternodes where they collide. The periaxonal width that stops internodal wave propagation (about 400 nm) is significantly larger than the highest value of the physiological range for this parameter (30 nm). Internodal activation is directly manifested as transmembrane internodal potential or as a full-sized action potential in periaxonal space where it can hardly be detected, and only as a small deflection in intracellular space. However, changes in the periaxonal potential cause transmyelin currents that lead to significant DAP. The shape and amplitude of DAP depends on myelin parameters and densities of internodal channels. Several technical parameters affect the results of calculations. Internodal spatial segmentation has to be sufficiently fine (at most 20 microm) for the model to be able to simulate internodal activation. We employ 338 internodal segments as compared with up to 21 used in previous models. Ionic accumulation together with related diffusive and electrical processes alter the calculated DAP amplitude. Inclusion of these processes in calculations demands such increase in the total number of segments that the numerical methods used up to now become unapplicable. To overcome the problem, an iterative implicit approach is proposed. It reduces a matrix of general type in multi-cable models to tridiagonal one and accelerates calculations considerably.

Analysis of real-time numerical integration methods applied to dynamic clamp experiments

Real-time systems are frequently used as an experimental tool, whereby simulated models interact in real time with neurophysiological experiments. The most demanding of these techniques is known as the dynamic clamp, where simulated ion channel conductances are artificially injected into a neuron via intracellular electrodes for measurement and stimulation. Methodologies for implementing the numerical integration of the gating variables in real time typically employ first-order numerical methods, either Euler or exponential Euler (EE). EE is often used for rapidly integrating ion channel gating variables. We find via simulation studies that for small time steps, both methods are comparable, but at larger time steps, EE performs worse than Euler. We derive error bounds for both methods, and find that the error can be characterized in terms of two ratios: time step over time constant, and voltage measurement error over the slope factor of the steady-state activation curve of the voltage-dependent gating variable. These ratios reliably bound the simulation error and yield results consistent with the simulation analysis. Our bounds quantitatively illustrate how measurement error restricts the accuracy that can be obtained by using smaller step sizes. Finally, we demonstrate that Euler can be computed with identical computational efficiency as EE.

Predicted profiles of ion concentrations in olfactory cilia in the steady state

The role of ciliary geometry for transduction events was explored by numerical simulation. The changes in intraciliary ion concentrations, suspected to occur during transduction, could thus be estimated. The case of a single excised cilium, having a uniform distribution of membrane channels, voltage clamped to -80 mV, was especially investigated. The axial profile of membrane voltage was that of a leaky cable. The Ca(2+) concentration profile tended to show a maximum in proximal segments, due to a preponderance of Ca(2+) inflow over Ca(2+) export at those locations. The local increase in Ca(2+) concentration activated Cl(-) channels. The resulting current caused a local drop in Cl(-) concentration, especially at the tip of the cilium and in distal segments, accompanied by a drop in ciliary K(+) concentration. In consequence, the membrane Cl(-) current was low in distal segments but stronger in proximal segments, where resupply was sufficient. The model predicts that the Cl(-) depletion will codetermine the time course of the receptor potential or current and the ciliary stimulus-response curve. In conclusion, when modeling with transduction elements presently known to participate, the ciliary geometry has large effects on ion distributions and transduction currents because ciliary ion transport is limited by axial electrodiffusion.