calcium waves¶

calcium waves

Overview

In every cell type, intracellular Ca²⁺ signaling plays a major role in a diverse set of activities, including process homeostasis and genetic expression. In neurons, Ca²⁺ wave propagation has been implicated as a major component of intracellular signaling and is complementary to membrane electrical signaling via action potentials. It follows that the cell must contain tightly-controlled intrinsic mechanisms for regulating intracellular Ca²⁺ content and waves; these include cytosolic buffers, storage in the ER, and sequestration into mitochondria. In this tutorial, we examine how some of these basic mechanisms interact to produce a Ca²⁺ wave.

Concept

Ca²⁺ is well-buffered by several cytosolic proteins (e.g. calmodulin, calbindin, etc...), but must still traverse a significant distance within the neuron to reach its targets. This distance is long enough that pure diffusion would pose a temporal problem for cell communication. Instead, the cell uses oscillating Ca²⁺ waves to harness the power of chemical signaling. These waves can travel up and down the Endoplasmic Reticulum (ER) to effect homeostatic changes. Understanding how these waves are modulated in real time by the cell poses a great difficulty because of complex spatiotemporal patterns and nonuniform distribution of the channels that control Ca²⁺ flux across the ER. Inositol triphosphate (IP₃) receptors IP₃Rs), Sarco/endoplasmic reticulum Ca²⁺-ATPase (SERCA pump), and Ryanoine receptors (RyRs). For simplicity’s sake, we will not deal with RyRs in this tutorial.

Diffusible Ca²⁺ is released directly from the ER via IP₃Rs upon cooperative binding of free Ca²⁺ and IP₃ to the receptor. When IP₃ levels are ideal, this results in Ca²⁺-induced Ca²⁺ release (CICR), yielding a wave that can travel along the ER. As mentioned above, these waves are implicated in the gene expression and protein upregulation associated with neuronal plasticity by way of their interaction with the neuronal soma. The SERCA pumps within ER are normally responsible for maintaining appropriate intracellular Ca²⁺ levels. Two models exist for wave propagation: continuous and saltatory. To understand the physical and functional differences between these models, it is ideal to compare them using Reaction-Diffusion (RxD) wtihin the NEURON environment.

The model presented in this tutorial generates Ca²⁺ waves and is based loosely off of the model we used in Neymotin et al., 2015.

Specification

Python is a rich and user-friendly programming language that offers significant options for systems integration. For this reason, we will run the NEURON and RxD modules within the Python interpreter.

Like other RxD problems within NEURON, the first step is to import and load the appropriate libraries into the Python interpreter.

In [1]:

from neuron import h, rxd

It is often convenient to load NEURON's stdrun library, which provides high level functions for controlling integration; we do so here:

In [2]:

h.load_file('stdrun.hoc')

Out[2]:

1.0

Basic structure

First, we define parameters for the section to be used. Like any imported module within Python, the imported libraries can be accessed by the h and rxd prefixes:

In [3]:

sec = h.Section(name='sec')
sec.L = 101
sec.diam = 1
sec.nseg = 101

The above code defines a section of a neuron with Length = 101 ({\mu}m), diameter = 1 ({\mu}m), and 101 discrete segments. By chunking the neuron into a large number of segments, we increase the resolution of the output at the expense of processing power. For a smaller number of segments, the opposite is true. In larger models, it is useful to represent nseg as a variable – this way, if changes need to be made during any number of simulations, the single variable can be altered and will affect all subsequent iterations of nseg. In this example, we only reference a single section, so there is no need to elaborate.

Enable variable time step

Next, we activate the variable time-step integration algorithm. This is a class of functions which are part of the SUNDIALS package and allow variable-order, variable time-step, and implicit solving of ordinary differential equations (ODE). This is a particularly useful tool for solving a system with multiple species, whose concentrations may change suddenly from stable conditions.

In [4]:

h.CVode().active(True)

Out[4]:

True

The argument of active specifies whether or not variable step integration will be used. To use fix step integration, pass the argument False instead.

Parameters

The diffusion coefficients caDiff and ip3Diff here denote the effective diffusion rate after buffering by native cytosol proteins is taken into account.

In [5]:

caDiff = 0.016
ip3Diff = 0.283
cac_init = 1.e-4
ip3_init = 0.1
gip3r = 12040
gserca = 0.3913
gleak = 6.020
kserca = 0.1
kip3 = 0.15
kact = 0.4
ip3rtau = 2000
fc = 0.83
fe = 0.17

Regions

An rxd.Region is the volume in which species exist and reactions can occur.

In [6]:

cyt = rxd.Region(h.allsec(),
                 nrn_region='i',
                 geometry=rxd.FractionalVolume(fc, surface_fraction=1),
                 name='cyt')

This uses the rxd.Region class to define the properties for the region on the immediate interior of the membrane (nrn_region = "i"). Electrophysiologically, the only regions that matter for this simulation are the spaces directly inside the membranes. The “i” marks that we are concerned with the inner region only. For any region, a geometry may be specified. For this model, we make the assumption that ER is roughly evenly-mixed throughout the neuron segment. That is, in any given chunk of neuron, we assume there is a constant fraction of that chunk that is ER and a constant fraction that is cytosol; this is indicated using an rxd.FractionalVolume. The geometry argument tells NEURON what to expect the region to actually look like (i.e. a cylinder has a specific 3D spatial representation and therefore requires specific diffusion modeling for ions within that particular geometry).

We set surface_fraction=1 because for the purposes of this model we suppose that at all points the plasma membrane is in contact with the cytosol (i.e. we assume the ER does not contact the plasma membrane and that all calcium influx from ion channels must first enter the cytosol before it can be sequestered into the ER).

name is an optional argument, but it can be queried later and some tools support using it to identify the region.

We define the ER similarly, except we omit the surface_fraction argument as it defaults to 0; i.e. that the ER does not contact the plasma membrane:

In [7]:

er = rxd.Region(h.allsec(), geometry=rxd.FractionalVolume(fe), name='er')

Finally, we define the membrane separating the cytosol from the ER. rxd.MultiCompartmentReactions crossing this membrane will be scaled by its surface area. Given that we are assuming that all chunks of neuron have ER and cytosol mixed together in the same way, it follows that the area of the membrane within a volume should be proportional to the size of that volume. We indicate this using the rxd.DistributedBoundary geometry. This is typically the right rule for separating two volumes defined as rxd.FractionalVolumes:

In [8]:

cyt_er_membrane = rxd.Region(h.allsec(),
                             geometry=rxd.DistributedBoundary(2),
                             name='cyt_er_membrane')

Here the argument 2 to rxd.DistributedBoundary indicates that there is 2 μm² of ER boundary per cubic micron of volume. Ideally, this number would be based on empirical data, but for now we choose it arbitrarily.

Species

Recall that the definition of an rxd.Species has the general form:

In [9]:

s = rxd.Species(regions=None,
               d=0,
               name=None,
               charge=0,
               initial=None,
               atolscale=1)

where:

regions can either be a single region or a list of regions.
d is the diffusion coefficient for the species; if none is specified, the species does not diffuse.
name assigns a name to the Species that syncs with the rest of NEURON (NMODL, GUI tools, variables of the form e.g. cai, etc)
charge refers to the charge of the Species; e.g. a calcium ion would have charge 2
initial refers to the initial concentration of the Species. If initial is not specified and the Species has a name and exists on an rxd.Region with nrn_region set, then NEURON's default initialization rules apply; otherwise the default is 0.
atolscale indicates a tolerance scale factor for variable step integration and is especially important for calcium and other substances with a low concentration

Because these are keyword arguments, they can be left out or even entered in any order.

We suppose calcium exists in both the ER and cytosol, and IP3 is present in the cytosol alone:

In [10]:

caCYT_init = 0.1
ca = rxd.Species([cyt, er], d=caDiff, name='ca', charge=2, initial=caCYT_init)
ip3 = rxd.Species(cyt, d=ip3Diff, initial=ip3_init)

We apologize, but this tutorial is incomplete and under revision. Please continue to explore the other examples. If you want to get a sense of where this is going, this is based on an older example but that file has a number of issues.