neuron.rxd.generalizedReaction.GeneralizedReaction

neuron.rxd.geometry.FixedCrossSection
__call__ · __init__

neuron.rxd.geometry.FixedPerimeter
__call__ · __init__

neuron.rxd.geometry.FractionalVolume
__call__ · __init__

neuron.rxd.geometry.Shell
__call__ · __init__

neuron.rxd.MultiCompartmentReaction

neuron.rxd.node.Node
__init__ · _ref_concentration · concentration · d · satisfies

neuron.rxd.nodelist.NodeList
__call__ · __init__ · concentration · diff · region · species · surface_area · volume · x

neuron.rxd.Rate
__init__

neuron.rxd.Reaction
__init__

neuron.rxd.Region
__init__ · nrn_region · secs

neuron.rxd.section1d.Section1D
__init__ · L · name · nrn_region · nseg · region · section_orientation · species

neuron.rxd.Species
__init__ · charge · indices · name · nodes

neuron.rxd.species.SpeciesOnRegion
__init__ · concentration · indices · nodes

Basic Reaction-Diffusion

Warning

These functions are experimental. Please contact the NEURON developers with issues and be sure to validate your results.

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

Declare a conceptual region of the neuron.

Examples: Cytosol, ER

__init__(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.

Note

dimension and dx will be deprecated in a future version

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

Note

Setting only supported in 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.

Note

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


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

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

Declare a species.

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 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)

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.

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

property charge

Get or set the charge of the Species.

Note

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

indices(r=None, secs=None)

return the indices corresponding to this species in the given region

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.

property name

Get or set the name of the Species.

This is used only for syncing with the non-reaction-diffusion parts of NEURON.

Note

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

property nodes

A NodeList of all the nodes corresponding to the species.

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


class neuron.rxd.Reaction(*args, **kwargs)
__init__(*args, **kwargs)

Specify a reaction to be added to the system.

Examples:
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

)


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

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

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.

__init__(species, rate, regions=None, membrane_flux=False)

create a rate of change for a species on a given region or set of regions

if regions is None, then does it on all regions


class neuron.rxd.MultiCompartmentReaction(*args, **kwargs)

class neuron.rxd.node.Node
__init__()

Initialize self. See help(type(self)) for accurate signature.

property _ref_concentration

Returns a HOC reference to the Node’s concentration

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

property concentration

Gets the concentration at the Node.

property d

Gets the diffusion rate within the compartment.

satisfies(condition)

Tests if a Node satisfies a given condition.

If a nrn.Section object or RxDSection is provided, returns True if the Node lies in the section; else False. If a Region object is provided, returns True if the Node lies in the Region; else False. If a number between 0 and 1 is provided, returns True if the normalized position lies within the Node; else False.


class neuron.rxd.generalizedReaction.GeneralizedReaction

an abstract class, parent of Rate, Reaction, MultiCompartmentReaction


class neuron.rxd.species.SpeciesOnRegion(species, region)
__init__(species, region)

The restriction of a Species to a Region.

property concentration

An iterable of the current concentrations.

indices(r=None, secs=None)

If no Region is specified or if r is the Region specified in the constructor, returns a list of the indices of state variables corresponding to the Species when restricted to the Region defined in the constructor.

If r is a different Region, then returns the empty list. If secs is a set return all indices on the regions for those sections.

property nodes

A NodeList of the Node objects containing concentration data for the given Species and Region.

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 Node objects belonging to the restriction ca[cyt] whereas the second option constructs all Node objects belonging to the Species ca and then culls the list to only include those also belonging to the Region cyt.


class neuron.rxd.section1d.Section1D(species, sec, diff, r)
property L

The length of the section in microns.

__init__(species, sec, diff, r)

Initialize self. See help(type(self)) for accurate signature.

name()

The name of the section, as defined by HOC.

property nrn_region

The HOC region, if any.

‘i’ if the corresponding Region maps to HOC’s internal concentration

‘o’ if it maps to HOC’s external concentration

None otherwise.

property nseg

The number of segments in the NEURON section.

In some implementations, this might be different than the number of nodes in the section.

property region

The Region containing the RxDSection.

property section_orientation

The HOC section orientation.

property species

The Species object stores on this RxDSection.


class neuron.rxd.nodelist.NodeList(items)
__call__(restriction)

returns a sub-NodeList consisting of nodes satisfying restriction

__init__(items)

Constructs a NodeList from items, a python iterable containing Node objects.

property concentration

Returns the concentration of the Node objects in the NodeList as an iterable.

property diff

Returns the diffusion constant of the Node objects in the NodeList as an iterable.

property region

An iterable of the Region objects corresponding to the Node objects in the NodeList.

Read only.

property species

An iterable of the Species objects corresponding to 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 volume

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

Read only.

property x

An iterable of the normalized positions of the Node objects in the NodeList.

Read only.


class neuron.rxd.geometry.FixedCrossSection(cross_area, surface_area=0)
__call__()

calling returns self to allow for rxd.inside or rxd.inside()

__init__(cross_area, surface_area=0)

Initialize self. See help(type(self)) for accurate signature.


class neuron.rxd.geometry.FractionalVolume(volume_fraction=1, surface_fraction=0, neighbor_areas_fraction=None)
__call__()

calling returns self to allow for rxd.inside or rxd.inside()

__init__(volume_fraction=1, surface_fraction=0, neighbor_areas_fraction=None)

Initialize self. See help(type(self)) for accurate signature.


class neuron.rxd.geometry.FixedPerimeter(perimeter, on_cell_surface=False)
__call__()

calling returns self to allow for rxd.inside or rxd.inside()

__init__(perimeter, on_cell_surface=False)

Initialize self. See help(type(self)) for accurate signature.


class neuron.rxd.geometry.Shell(lo=None, hi=None)
__call__()

calling returns self to allow for rxd.inside or rxd.inside()

__init__(lo=None, hi=None)

Initialize self. See help(type(self)) for accurate signature.