rxd.export.sbml · rxd.nthread · rxd.re_init · rxd.set_solve_type

rxd.DistributedBoundary

rxd.Extracellular

rxd.FixedCrossSection

rxd.FixedPerimeter

rxd.FractionalVolume

rxd.MultiCompartmentReaction

rxd.node.Node
include_flux · satisfies

rxd.node.Node1D
satisfies

rxd.node.Node3D
satisfies

rxd.node.NodeExtracellular
satisfies

rxd.nodelist.NodeList
include_flux · value_to_grid

rxd.Parameter

rxd.Rate

rxd.Reaction

rxd.Region

rxd.RxDException

rxd.ScalableBorder

rxd.Shell

rxd.Species
concentrations · defined_on_region · indices · re_init

rxd.State

Basic Reaction-Diffusion

Overview

NEURON provides the rxd submodule to simplify and standardize the specification of models incorporating reaction-diffusion dynamics, including ion accumulation. To load this module, use:

from neuron import rxd

Note: In older code, you may find from neuron import crxd as rxd but this is equivalent to the above as crxd has been an alias for rxd for several years.

In general, a reaction-diffusion model specification involves answering three conceptual questions:

  1. Where the dynamics are occurring (specified using an rxd.Region or rxd.Extracellular)
  2. Who is involved (specified using an rxd.Species or rxd.State)
  3. What the reactions are (specified using rxd.Reaction, rxd.Rate, or rxd.MultiCompartmentReaction)

Another key class is rxd.Parameter for defining spatially varying parameters. Integration options may be specified using rxd.set_solve_type().

Specifying the domain

NEURON provides two main classes for defining the domain where a given set of reaction-diffusion rules applies: rxd.Region and rxd.Extracellular for intra- and extracellular domains, respectively. Once defined, they are generally interchangeable in the specification of the species involved, the reactions, etc. The exact shape of intracellular regions may be specified using any of a number of geometries, but the default is to include the entire intracellular space.

Intracellular regions and regions in Frankenhauser-Hodgkin space

class rxd.Region

Declares a conceptual intracellular region of the neuron.

Syntax:

r = rxd.Region(secs=None, nrn_region=None, geometry=None, dimension=None, dx=None, name=None)

In NEURON 7.4+, secs is optional at initial region declaration, but it must be specified before the reaction-diffusion model is instantiated.

Here:

  • secs is any Python iterable of sections (e.g. soma.wholetree() or [soma, apical, basal] or h.allsec())
  • nrn_region specifies the classic NEURON region associated with this object and must be either "i" for the region just inside the plasma membrane, "o" for the region just outside the plasma membrane or None for none of the above.
  • name is the name of the region (e.g. cyt or er); this has no effect on the simulation results but it is helpful for debugging
  • dx deprecated; do not use
  • dimension deprecated; do not use
property nrn_region

Get or set the classic NEURON region associated with this object.

There are three possible values:

  • 'i' – just inside the plasma membrane
  • 'o' – just outside the plasma membrane
  • None – none of the above

Setting requires NEURON 7.4+, and then only before the reaction-diffusion model is instantiated.

property secs

Get or set the sections associated with this region.

The sections may be expressed as a NEURON SectionList or as any Python iterable of sections.

Note: The return value is a copy of the internal section list; modifying
it will not change the Region.

Setting requires NEURON 7.4+ and allowed only before the reaction-diffusion model is instantiated.

property geometry

Get or set the geometry associated with this region.

Setting the geometry to None will cause it to default to rxd.geometry.inside.

Added in NEURON 7.4. Setting allowed only before the reaction-diffusion model is instantiated.

property name

Get or set the Region’s name.

Added in NEURON 7.4.

For specifying the geometry

NEURON provides several built-in geometries for intracellular simulation that may be specified by passing a geometry argument to the rxd.Region constructor. New region shapes may be defined by deriving from neuron.rxd.geometry.RxDGeometry.

rxd.inside

The entire region inside the cytosol. This is the default. Use via e.g.

cyt = rxd.Region(h.allsec(), name="cyt", nrn_region="i", geometry=rxd.inside)
rxd.membrane

The entire plasma membrane.

cyt = rxd.Region(h.allsec(), name="cyt", nrn_region="i", geometry=rxd.membrane)
class rxd.DistributedBoundary

Boundary that scales with area.

b = rxd.DistributedBoundary(area_per_vol, perim_per_area=0)

area_per_vol is the area of the boundary as a function of the volume containing it. e.g.

g = rxd.DistributedBoundary(2) defines a geometry with 2 um^2 of area per every um^3 of volume.

perim_per_area is the perimeter (in µm) per 1 µm^2 cross section of the volume. For use in reaction-diffusion problems, it may be safely omitted if and only if no species in the corresponding region diffuses.

This is often useful for separating rxd.FractionalVolume objects.

It is assumed that the area is always strictly on the interior.

class rxd.FractionalVolume

Defines a geometry occupying a set fraction of the cross-sectional area. e.g. perhaps the cytosol would occupy 0.83 of the cross-section (and all the surface) but the ER would only occupy 0.17 of the cross-section and none of the surface.

Syntax:

r = rxd.FractionalVolume(volume_fraction=1, surface_fraction=0, neighbor_areas_fraction=None)
class rxd.Shell

Defines a radial shell inside or outside of the plasma membrane. This is sometimes used to simulate a 2D-style diffusion where molecules move both longitudinally and into/out of the center of the dendrite. Consider using 3D simulation instead, which is better able to represent branch point dynamics.

Syntax:

r = rxd.Shell(lo=None, hi=None)

Example:

See the radial diffusion tutorial.
class rxd.FixedCrossSection

Syntax:

r = rxd.FixedCrossSection(cross_area, surface_area=0)
class rxd.FixedPerimeter

Syntax:

r = rxd.FixedPerimeter(perimeter, on_cell_surface=False)
class rxd.ScalableBorder

A membrane that scales proportionally with the diameter

Example use:

  • the boundary between radial shells

Sometimes useful for the boundary between rxd.FractionalVolume objects, but see also rxd.DistributedBoundary which scales with area.

Extracellular regions

class rxd.Extracellular

Declare a extracellular region to be simulated in 3D; unlike rxd.Region, this always describes the extracellular region.

Syntax:

r = rxd.Extracellular(xlo, ylo, zlo, xhi, yhi, zhi, dx,
                      volume_fraction=1, tortuosity=None, permeability=None)

Here:

  • xlo, …, zhi – define the bounding box of the domain; it should typically contain the entire morphology… by default NEURON assumes reflective (Neumann) boundary conditions, so you may want the box to extend well beyond the cell morphology depending on your use case
  • dx – voxel edge size in µm
  • tortuosity – increase factor in path length due to obstacles, effective diffusion coefficient d/tortuosity^2; either a single value for the whole region or a Vector giving a value for each voxel. Default is 1 (no change).
  • volume_fraction – the free fraction of extracellular space; a volume_fraction of 1 assumes no cells; lower values are probably warranted for most simulations

Example:

A tutorial demonstrating extracellular diffusion is available here.

Defining proteins, ions, etc

Values that are distributed spatially on an rxd.Region or rxd.Extracellular may be defined using an rxd.Species if they represent things that change and diffuse, an rxd.State if they’re in fixed locations but changeable (e.g. gates in an IP3R), or an rxd.Parameter if they are just fixed values.

class rxd.Species

Declare an ion/protein/etc that can react and diffuse.

Syntax:

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

Parameters:

  • regions – a Region or list of Region objects containing the species
  • d – the diffusion constant of the species (optional; default is 0, i.e. non-diffusing)
  • name – the name of the Species; used for syncing with NMODL and HOC (optional; default is none)
  • charge – the charge of the Species (optional; default is 0)
  • initial – the initial concentration or None (if None, then imports from HOC if the species is defined at finitialize, else 0); can be specified as a constant or as a function of an rxd.Node
  • atolscale – scale factor for absolute tolerance in variable step integrations
  • ecs_boundary_conditions – if Extracellular rxd is used ecs_boundary_conditions=None for zero flux boundaries or if ecs_boundary_conditions=the concentration at the boundary.
  • represents – optionally provide CURIE (Compact URI) to annotate what the species represents e.g. CHEBI:29101 for sodium(1+)

Note

Charge must match the charges specified in NMODL files for the same ion, if any. Common species charges include: sodium (+1), potassium (+1), calcium (+2), magnesium (+2), chloride (-1).

You probably want to adjust atolscale for species present at low concentrations (e.g. calcium).

NEURON does not require any specific ontology for identifiers, however CHEBI contains identifiers for many substances of interest in reaction-diffusion modeling. A number of ontology search providers are available on the internet, including BioPortal.

To refer to a given Species restricted to a certain region, specify the region in square brackets. e.g. er_calcium = ca[er].

property nodes

An rxd.nodelist.NodeList of all the nodes corresponding to the species.

This can then be further restricted using the callable property of NodeList objects.

Example: All nodes from sp on the Section dend:

nodes_on_dend = sp.nodes(dend)

Example: All nodes from sp on the Segment dend(0.5):

nodes_on_dend_center = sp.nodes(dend(0.5))

For 1D simulation with a species defined on 1 region, this will be a rxd.nodelist.NodeList of length 1 (or 0 if the species is not defined on the segment); for 3D simulation or with multiple regions, the list may be longer and further filtering may be required.

Note

The code

node_list = ca[cyt].nodes

is more efficient than the otherwise equivalent

node_list = ca.nodes(cyt)

because the former only creates the rxd.Node objects belonging to the restriction ca[cyt] whereas the second option constructs all rxd.Node objects belonging to the rxd.Species ca and then culls the list to only include those also belonging to the rxd.Region cyt.

property states

A list of a copy of the state values corresponding to this species; modifying the values in the list will not change the values used in the simulation.

re_init()

Syntax:

sp.re_init()

where sp is an instance of an rxd.Species.

Reinitialize the rxd concentration of this species to match the NEURON grid. Used when e.g. cai is modified directly instead of through a corresponding rxd.Species.

concentrations()

Deprecated. Do not use.

property charge

Get or set the charge of the Species.

Setting was added in NEURON 7.4 and is allowed only before the reaction-diffusion model is instantiated.

property regions

Get or set the regions where the Species is present.

Setting was added in NEURON 7.4 and is allowed only before the reaction-diffusion model is instantiated. Getting was added in NEURON 7.5.

indices()

Return the indices corresponding to this species in the given region.

This is occasionally useful, but it’s generally best to avoid explicit indexing when developing models.

Syntax:

sp.indices(r=None, secs=None)

where sp is an instance of an rxd.Species.

If r is None, then returns all species indices. If r is a list of regions return indices for only those sections that are on all the regions. If secs is a set return all indices on the regions for those sections.

defined_on_region()

Returns True if the rxd.Species sp is present on r, else False.

Syntax:

result = sp.defined_on_region(r)
property name

Get or set the name of the Species.

This is used only for syncing with the non-reaction-diffusion parts of NEURON (NMODL, HOC).

Setting requires NEURON 7.4+, and then only before the reaction-diffusion model is instantiated.

class rxd.State

Like an rxd.Species but indicates the semantics of something that is not intended to diffuse.

class rxd.Parameter

Declares a parameter, that can be used in place of a rxd.Species, but unlike rxd.Species a parameter will not change.

Syntax:

s = rxd.Parameter(regions, name=None, charge=0,
                  value=None, represents=None)

Parameters:

  • regions – a rxd.Region or list of rxd.Region objects containing the parameter
  • name – the name of the parameter; used for syncing with NMODL and HOC (optional; default is none)
  • charge – the charge of the Parameter (optional, probably only rarely needed; default is 0)
  • value – the value or None (if None, then imports from HOC if the parameter is defined at finitialize(), else 0)
  • represents – optionally provide CURIE (Compact URI) to annotate what the parameter represents e.g. CHEBI:29101 for sodium(1+)

Note

Charge must match the charges specified in NMODL files for the same ion, if any.

Attempting to specify a non-zero diffusion rate for an rxd.Parameter will raise an rxd.RxDException.

Defining reactions

NEURON provides three classes for specifying reaction dynamics: rxd.Reaction for single-compartment (local) reactions; rxd.MultiCompartmentReaction for reactions spanning multiple compartments (e.g. a pump that moves calcium from the cytosol into the ER changes concentration in two regions), and rxd.Rate for specifying changes to a state variable directly by an expression to be added to a differential equation. Developers may introduce new forms of reaction specification by subclassing neuron.rxd.generalizedReaction.GeneralizedReaction, but this is not necessary for typical modeling usage.

It is sometimes necessary to build rate expressions including mathematical functions. To do so, use the functions defined in neuron.rxd.rxdmath as listed below.

class rxd.Reaction

Syntax:

r1 = rxd.Reaction(reactant_sum, product_sum, forward_rate,
                  backward_rate=0, regions=region_list, custom_dynamics=False)

Here:

  • reactant_sum is a combination of rxd.Species, rxd.State, rxd.Parameter, e.g. ca + 2 * cl representing the left-hand-side of the reaction.
  • product_sum is like reactant_sum but represdenting the right-hand-side.
  • forward_rate and backward_rate represent the reaction rates; reactions are assumed to be governed by mass-action kinetics with these as the rate constants unless custom_dynamics is true, in which case these are expressions fully defining the rate of change. In particular, these can be constants or expressions combining rxd.Species etc with constants, arithmetic, and rxd.v
  • region_list is a list of regions on which this reaction occurs. If ommitted or None, the reaction occurs on all regions where all involved species are defined.

Examples:

For the following, we assume that H, O, H2O, dimer, and decomposed are instances of rxd.Species.

For 2 * H + O > H2O in a mass action reaction at rate k:

r = rxd.Reaction(2 * H + O, H2O, k)

To constrain the reaction to a specified list of regions, say to just the extracellular space (ext) and the cytosol (cyt), use the regions keyword, e.g.

r = rxd.Reaction(2 * H + O, H2O, k, regions=[ext, cyt])

For a bi-directional reaction, specify a backward reaction rate. e.g. if kf is the forward rate and kb is the backward rate, then:

r = rxd.Reaction(2 * H + O, H2O, kf, kb)

To use dynamics other than mass-action, add that mass_action=False flag and put the full formula instead of a mass-action rate for kf (and kb). E.g. Michaelis-Menten degradation

r = rxd.Reaction(
    dimer, decomposed, dimer / (k + diamer), mass_action=False
)
property f_rate

Get or set the forward reaction rate.

property b_rate

Get or set the backward reaction rate.

class rxd.Rate

Declare a contribution to the rate of change of a species or other state variable.

Syntax:

r = rxd.Rate(species, rate, regions=None, membrane_flux=False)

If regions is None, then the rate applies on all regions.

Example:

constant_production = rxd.Rate(protein, k)

If this was the only contribution to protein dynamics and there was no diffusion, the above would be equivalent to:

dprotein/dt = k

If there are multiple rxd.Rate objects (or an rxd.Reaction, etc) acting on the same species, then their effects are summed.

class rxd.MultiCompartmentReaction

Specify a reaction spanning multiple regions to be added to the system. Use this for, for example, pumps and channels, or interactions between species living in a volume (e.g. the cytosol) and species on a membrane (e.g. the plasma membrane). For each species/state/parameter, you must specify what region you are referring to, as it could be present in multiple regions. You must also specify a membrane or a border (these are treated as synonyms) that separates the regions involved in your reaction. This is necessary because the default behavior is to scale the reaction rate by the border area, as would be expected if one of the species involved is a pump that is binding to a species in the volume. If this is not the desired behavior, pass the keyword argument scale_by_area=False. Pass in membrane_flux=True if the reaction produces a current across the plasma membrane that should affect the membrane potential. Unlike rxd.Reaction objects, the base units for the rates are in terms of molecules per square micron per ms.

Mathematical functions for rate expressions

NEURON’s neuron.rxd.rxdmath module provides a number of mathematical functions that can be used to define reaction rates. These generally mirror the functions available through Python’s math module but support rxd.Species objects.

To use any of these, first do:

from neuron.rxd import rxdmath

Example:

degradation_switch = (1 + rxdmath.tanh((ip3 - threshold) * 2 * m)) / 2
degradation = rxd.Rate(ip3, -k * ip3 * degradation_switch)

For a full runnable example, see this tutorial which as here uses the hyperbolic tangent as an approximate on/off switch for the reaction.

Full documentation on this submodule is available here or you may go directly to the documentation for any of its specific functions: rxdmath.acos(), rxdmath.acosh(), rxdmath.asin(), rxdmath.asinh(), rxdmath.atan(), rxdmath.atan2(), rxdmath.ceil(), rxdmath.copysign(), rxdmath.cos(), rxdmath.cosh(), rxdmath.degrees(), rxdmath.erf(), rxdmath.erfc(), rxdmath.exp(), rxdmath.expm1(), rxdmath.fabs(), rxdmath.factorial(), rxdmath.floor(), rxdmath.fmod(), rxdmath.gamma(), rxdmath.lgamma(), rxdmath.log(), rxdmath.log10(), rxdmath.log1p(), rxdmath.pow(), rxdmath.pow(), rxdmath.sin(), rxdmath.sinh(), rxdmath.sqrt(), rxdmath.tan(), rxdmath.tanh(), rxdmath.trunc(), rxdmath.vtrap()

Manipulating nodes

A rxd.node.Node represents a particular state value or rxd.Parameter in a particular location. Individual rxd.node.Node objects are typically obtained either from being passed to an initialization function or by filtering or selecting from an rxd.nodelist.NodeList returned by rxd.Species.nodes. Node objects are often used for recording concentration using rxd.node.Node._ref_concentration.

class rxd.nodelist.NodeList

An rxd.nodelist.NodeList is a subclass of a Python list containing rxd.node.Node objects. It is not intended to be created directly in a model, but rather is returned by rxd.Species.nodes.

Standard Python list methods are supported, including .append(node), .extend(node_list), .insert(i, node), .index(node), and manipulation of lists like len(node_list), node_list[0], or node_list[5:12]. Additionally one may iterate over a NodeList as in:

for node in ca.nodes:
    ...

(Here ca is assumed to be an rxd.Species and thus ca.nodes is an rxd.nodelist.NodeList.)

A key added functionality is the ability to filter the nodes by rxd property (returning a new rxd.nodelist.NodeList). Any filter object supported by the .satifies method of the node types present in the rxd.nodelist.NodeList may be passed in parentheses; e.g.

To filter the rxd.Species ca’s nodes for just the ones present on the nrn.Segment dend(0.5), use:

new_node_list = ca.nodes(dend(0.5))

To filter the new_node_list to only contain nodes present in the rxd.Region er:

just_er = new_node_list(er)

In addition, the following methods and properties are supported:

property value

Gets or sets the values associated with the stored nodes. Getting always returns a list, even if the rxd.nodelist.NodeList has length 0 or 1. Setting may be to a constant (in which case all nodes are set to the same value) or to a list (in which case the list values are assigned in order to the nodes). In the latter case, if the length of the list does not match the length of the node list, an rxd.RxDException is raised.

The list that is returned by reading this property is a copy of the underlying data; that is, changing it will have no effect on the values stored.

This currently has the same behavior as rxd.nodelist.NodeList.concentration however in the future these are intended to be different for stochastic simulation.

property concentration

Gets or sets the concentration associated with the stored nodes. Getting always returns a list, even if the rxd.nodelist.NodeList has length 0 or 1. Setting may be to a constant (in which case all nodes are set to the same value) or to a list (in which case the list values are assigned in order to the nodes). In the latter case, if the length of the list does not match the length of the node list, an rxd.RxDException is raised.

The list that is returned by reading this property is a copy of the underlying data; that is, changing it will have no effect on the values stored.

This currently has the same behavior as rxd.nodelist.NodeList.value however in the future these are intended to be different for stochastic simulation.

property segment

Returns a list of the nrn.Segment objects associated with the nodes in the NodeList.

The list that is returned by reading this property is a copy of the underlying data; that is, changing it will have no effect on the values stored.

property _ref_value

A pointer to the memory location storing the rxd.node.Node.value when the NodeList has length 1; otherwise an rxd.RxDException is raised.

property _ref_concentration

A pointer to the memory location storing the rxd.node.Node.concentration when the NodeList has length 1; otherwise an rxd.RxDException is raised.

property diff

Get or set the diffusion constants of the contained Node objects. Getting returns a list that is a copy of the underlying data. Setting accepts either a constant or a list of matching length; passing a list of a different length raises an rxd.RxDException.

property volume

An iterable of the volumes of the Node objects in the NodeList.

Read only.

property surface_area

An iterable of the surface areas of the Node objects in the NodeList.

Read only.

property region

An iterable of the rxd.Region (or rxd.Extracellular) objects of the Node objects in the NodeList.

Read only.

property species

An iterable of the rxd.Species (or rxd.State or rxd.Parameter, as appropriate) objects of the Node objects in the NodeList.

Read only.

property x

An iterable of the normalized positions of the Node objects in the NodeList. Note: these values are always between 0 and 1 and represent positions within the corresponding Section. For 3D position, query the x3d property of the rxd.node.Node objects themselves.

Read only.

include_flux()

Includes the specified flux on all nodes in the NodeList. All arguments are passed directly to the underlying rxd.node.Node objects.

value_to_grid()

Returns a regular grid with the values of the 3D nodes in the list. This is sometimes useful for volumetric visualization however the generated array may be large in certain models.

The grid is a copy only.

Grid points not belonging to the object are assigned a value of NaN.

Nodes that are not 3d will be ignored. If there are no 3D nodes, returns a 0x0x0 numpy array.

Warning: Currently only supports nodelists over 1 region.

class rxd.node.Node

The combination of a single rxd.Species etc and a unique spatial location at whatever resolution (i.e. could be a segment and a region, or could be a 3D voxel and a region).

These objects are passed to an initialization function for rxd Species, States, and Parameters as ways of identifying a location. They are also useful for specifying localized fluxes or to record state variables.

There are three subtypes: rxd.node.Node1D, rxd.node.Node3D, and rxd.node.NodeExtracellular. They all support the methods and properties described here as well as some unique to their case features.

include_flux()

Include a flux contribution to a specific node.

The flux can be described as a NEURON reference, a point process and a property, a Python function, or something that evaluates to a constant Python float.

Supported units: molecule/ms, mol/ms, mmol/ms == millimol/ms == mol/s

Examples:

node.include_flux(mglur, 'ip3flux')           # default units: molecule/ms
node.include_flux(mglur, 'ip3flux', units='mol/ms') # units: moles/ms
node.include_flux(mglur._ref_ip3flux, units='molecule/ms')
node.include_flux(lambda: mglur.ip3flux)
node.include_flux(lambda: math.sin(h.t))
node.include_flux(47)

Warning:

Flux denotes a change in mass not a change in concentration. For example, a metabotropic synapse produces a certain amount of substance when activated. The corresponding effect on the node’s concentration depends on the volume of the node. (This scaling is handled automatically by NEURON’s rxd module.)
satisfies()

Tests if a Node satisfies a given condition.

Syntax:

result = node.satisfies(condition)

If a Section object or RxDSection is provided, returns True if the Node lies in the section; else False.

If a rxd.Region object is provided, returns True if the Node lies in the Region; else False.

Additional options are supported by specific subclasses, see rxd.node.Node1D.satisfies(), rxd.node.Node3D.satisfies(), and rxd.node.NodeExtracellular.satisfies().

property _ref_concentration

Returns a NEURON reference to the Node’s concentration. This result is typically passed to Vector.record() to record the concentration changes at a location over time.

(The node must store concentration data. Use rxd.node.Node._ref_molecules for nodes storing molecule counts.)

property _ref_molecules

Returns a NEURON reference to the Node’s concentration

(The node must store molecule counts. Use _ref_concentrations for nodes storing concentration.)

property _ref_value

Returns a NEURON reference to the Node’s value. This method always works, regardless of if the node stores a concentration or not.

property d

Get or set the diffusion rate within the compartment.

property concentration

Get or set the concentration at the Node.

Currently does not support nodes storing molecule counts. Use rxd.node.Node.molecules instead; attempting to use with a molecule count node will raise an rxd.RxDException.

property molecules

Get or set the number of molecules at the Node.

Currently does not support nodes storing concentrations. Use rxd.node.Node.concentration instead; attempting to use with a concentration node will raise an rxd.RxDException.

property value

Get or set the value associated with this Node.

For Species nodes belonging to a deterministic simulation, this is a concentration. For Species nodes belonging to a stochastic simulation, this is the molecule count.

property x3d

The 3D x-coordinate of the center of this Node.

property y3d

The 3D y-coordinate of the center of this Node.

property z3d

The 3D z-coordinate of the center of this Node.

property region

The rxd.Region or rxd.Extracellular containing the compartment. Read only.

property species

The rxd.Species, rxd.State, or rxd.Parameter containing the compartment. Read only.

property volume

The volume of the region spanned by the Node.

class rxd.node.Node1D

A subclass of rxd.node.Node used only for nodes being simulated in 1D.

satisfies()

Supports the options of rxd.node.Node.satisfies() and:

If a number between 0 and 1 is provided, returns True if the normalized position lies within the Node; else False.

property sec

The section containing the node. Read-only.

property segment

The segment containing the node. Read-only.

property x

The normalized position of the center of the node. Read-only.

property surface_area

The surface area of the compartment in square microns.

This is the area (if any) of the compartment that lies on the plasma membrane and therefore is the area used to determine the contribution of currents (e.g. ina) from mod files or KSChan to the compartment’s concentration.

Read only.

class rxd.node.Node3D

A subclass of rxd.node.Node used only for intracellular nodes being simulated in 3D.

satisfies()

Supports the options of rxd.node.Node.satisfies() and:

If a tuple is provided of length 3, return True if the Node contains the (x, y, z) point; else False.

property sec

The section containing the node. Read-only.

property segment

The segment containing the node. Read-only.

property surface_area

The surface area of the compartment in square microns.

This is the area (if any) of the compartment that lies on the plasma membrane and therefore is the area used to determine the contribution of currents (e.g. ina) from mod files or KSChan to the compartment’s concentration.

Read only.

class rxd.node.NodeExtracellular

A subclass of rxd.node.Node used only for extracellular nodes being simulated in 3D.

satisfies()

Supports the options of rxd.node.Node.satisfies() and:

If a tuple is provided of length 3, return True if the Node contains the (x, y, z) point; else False.

Membrane potential

property rxd.v

A special object representing the local membrane potential in a reaction-rate expression. This can be used with rxd.Rate and rxd.MultiCompartmentReaction to build ion channel models as an alternative to using NMODL, NeuroML (and converting to NMODL via jneuroml), the ChannelBuilder, or KSChan.

(If you want a numeric value for the current membrane potential at a segment seg use seg.v instead.)

Example (adapted from the Hodgkin Huxley via rxd tutorial)

from neuron.rxd.rxdmath import vtrap, exp, log
from neuron import rxd

alpha = 0.01 * vtrap(-(rxd.v + 55.0), 10.0)
beta = 0.125 * exp(-(rxd.v + 65.0)/80.0)
ntau = 1.0/(alpha + beta)
ninf = alpha/(alpha + beta)

# ... define cyt, mem, sections ...

ngate = rxd.State([cyt, mem], name='ngate', initial=0.24458654944007166)
n_gate = rxd.Rate(ngate, (ninf - ngate)/ntau)

Synchronization with segments

Changes to rxd.Species node concentrations are propagated to segment-level concentrations automatically no later than the next time step. This is generally the right direction for information to flow, however NEURON also provides a rxd.re_init() function to transfer data from segments to rxd.Species.

rxd.re_init()

Reinitialize all rxd.Species, rxd.State, and rxd.Parameter from changes made to NEURON segment-level concentrations. This calls the corresponding rxd.Species.re_init() methods. Note that reaction-diffusion models may contain concentration data at a finer-resolution than that of a nrn.Segment (e.g. for models being simulated in 3D).

Syntax:

rxd.re_init()

Numerical options

rxd.nthread()

Specify a number of threads to use for extracellular and 3D intracellular simulation. Currently has no effect on 1D reaction-diffusion models.

Syntax:

rxd.nthread(num_threads)

Example:

To simulate using 4 threads:

rxd.nthread(4)

Thread scaling performance is discussed in the NEURON extracellular and 3D intracellular methods papers.

rxd.set_solve_type()

Specify numerical discretization and solver options. Currently the main use is to indicate Sections where reaction-diffusion should be simulated in 3D.

Syntax:

rxd.set_solve_type(domain=None, dimension=None, dx=None, nsubseg=None, method=None)

where:

  • domain – a Section or Python iterable of sections. If the domain is None or omitted, the specification will apply to the entire model.
  • dimension – 1 or 3
  • dx – not implemented; specify dx for 3D models when creating the rxd.Region
  • nsubseg – not implemented
  • method – not implemented

This function may be called multiple times; the last setting for any given field will be used. Different sections may be simulated in different dimensions (a so-called hybrid model).

Warning

For 3D reaction-diffusion simulations, we recommend upgrading to at least NEURON 8.1.

(Calculation of 3D concentration changes from MOD file activity in NEURON 7.8.x and 8.0.x was underestimated due to an inconsistency in surface voxel partial volume calculations.)

SBML Export

rxd.export.sbml()

Export dynamics at a segment to an SBML file.

Syntax:

rxd.export.sbml(segment, filename=None, model_name=None, pretty=True)

This does not currently support rxd.MultiCompartmentReaction; attempting to export dynamics that involve such reactions will raise an rxd.RxDException.

Note

rxd.export is not available simply via from neuron import rxd; you must also:

import neuron.rxd.export

Saving and restoring state

Some simulations require a lengthy initialization before exploring various possible stimuli. In these situations, it is often convenient to run the initialization once, save the state, do an experiment, and revert back to the saved state.

Beginning in NEURON 8.1, reaction-diffusion states are included automatically when using SaveState which additionally saves many other model states.

If one wants to save and restore only reaction-diffusion states, this can be done via the following functions:

The reaction-diffusion state data returned by rxd.save_state() and expected by rxd.restore_state() consists of 16 bytes of metadata (8 bytes for a version identifier and 8 bytes for the length of the remaining portion) followed by gzip-compresssed state values. In particular, not every binary string of a given length is a valid state vector, nor is every state vector for a given model necessarily the same length (as the compressability may be different).

Error handling

Most errors in the usage of the neuron.rxd module should raise a rxd.RxDException.

class rxd.RxDException

An exception originating from the neuron.rxd module due to invalid usage. This allows distinguishing such exceptions from other errors.

The text message of an rxd.RxDException e may be read as str(e).