diffusion intracellular->extracellular->intracellular

[bremen]

Hi.

After reading the 2018 paper, i was wondering if is feasible, with rxd.extracellular, to reproduce the following situation:

A substance X is internally produced by cell_1.
It is pumped outside the membrane into the extracellular space.
Here it diffuses until it reaches cell_2.
It is then uptake by cell_2 and stored into an internal region.

Best regards
Stefano

[bremen]

Hi.
I was able to solve most of the workflow but I'm missing something in the "MultiCompartmentReaction".

I'm trying to release the substance X, correctly generated in the internal region of cell_1, to the ecs.
With the MultiCompartmentReaction wrote in this way:

Code: Select all

rxd.MultiCompartmentReaction(sub_x[int_region_cell1] > sub_x[ecs], 1, membrane = cell_1_membrane)

NEURON 7.7 (compiled 8 days ago) generates this error:

Code: Select all

Traceback (most recent call last):
  File "/usr/local/lib/python2.7/dist-packages/neuron/crxd/rxd.py", line 1240, in _init
    _compile_reactions()
  File "/usr/local/lib/python2.7/dist-packages/neuron/crxd/rxd.py", line 1074, in _compile_reactions
    creg._initalize()
  File "/usr/local/lib/python2.7/dist-packages/neuron/crxd/region.py", line 157, in _initalize
    self._ecs_initalize()
  File "/usr/local/lib/python2.7/dist-packages/neuron/crxd/region.py", line 123, in _ecs_initalize
    self.ecs_location_index[rid][sid][seg_idx] = s.index_from_xyz(x,y,z)
AttributeError: 'weakref' object has no attribute 'index_from_xyz'

Stefano

[adamjhn]

Thanks for reporting this bug, a fix has been pushed to the development branch. Currently the easiest way to do this is with mod files describing pump mechanism. The production, diffusion and accumulation are then handled by rxd. For example;

Code: Select all

from neuron import h, crxd as rxd
from matplotlib import pyplot

h.load_file('stdrun.hoc')
rxd.options.enable.extracellular = True

# create cell1 where `x` will be created and leak out
cell1 = h.Section('cell1')
cell1.pt3dclear()
cell1.pt3dadd(-20,0,0,10)
cell1.pt3dadd(-10,0,0,10)
cell1.insert('pump')

# create cell2 where `x` will be pumped in and accumulate
cell2 = h.Section('cell2')
cell2.pt3dclear()
cell2.pt3dadd(10,0,0,10)
cell2.pt3dadd(20,0,0,10)
cell2.insert('pump')

# Where?
# the intracellular spaces
cyt = rxd.Region(h.allsec(), name='cyt', nrn_region='i', geometry=rxd.FractionalVolume(0.9, surface_fraction=1.0))

org = rxd.Region(h.allsec(), name='org', geometry=rxd.FractionalVolume(0.1))

cyt_org_membrane = rxd.Region(h.allsec(), name='mem', geometry = rxd.ScalableBorder(1, on_cell_surface=False))


# the extracellular space
ecs = rxd.Extracellular(-55, -55, -55, 55, 55, 55, dx=10, volume_fraction=0.2, tortuosity=1.6)

# Who?
x = rxd.Species([cyt, org, cyt_org_membrane, ecs], name='x', d=1.0, charge=1, initial=0)
Xcyt = x[cyt]
Xorg = x[org]


# What? - produce X in cell 1
# parameter to limit production to cell 1
cell1_param = rxd.Parameter(org, initial=lambda node: 1.0 if node.segment.sec == cell1 else 0)

# production with a rate following Michaels Menton kinetics
createX = rxd.Rate(Xorg, cell1_param[org] * 1.0/(10.0 + Xorg))

# leak between organelles and cytosol
cyt_org_leak = rxd.MultiCompartmentReaction(Xcyt, Xorg, 1e4, 1e4, membrane=cyt_org_membrane)


# record the concentrations in the cells
t_vec = h.Vector()
t_vec.record(h._ref_t)
cell1_X = h.Vector()
cell1_X.record(cell1(0.5)._ref_xi)
cell1_Xorg = h.Vector()
cell1_Xorg.record(Xorg.nodes(cell1)(0.5)[0]._ref_concentration)

cell2_X = h.Vector()
cell2_X.record(cell2(0.5)._ref_xi)
cell2_Xorg = h.Vector()
cell2_Xorg.record(Xorg.nodes(cell2)(0.5)[0]._ref_concentration)


# run and plot the results
h.finitialize(-65)
h.continuerun(1000)

fig = pyplot.figure()
pyplot.subplot(2, 2, 1)
pyplot.plot(t_vec,cell1_X, label='cyt')
pyplot.plot(t_vec,cell1_Xorg, label='org')
pyplot.legend()
pyplot.xlabel('t (ms)')
pyplot.ylabel('cell 1 x (mM)')

pyplot.subplot(2, 2, 2)
pyplot.plot(t_vec,cell2_X, label='cyt')
pyplot.plot(t_vec,cell2_Xorg, label='org')
pyplot.xlabel('t (ms)')
pyplot.ylabel('cell 2 x (mM)')

pyplot.subplot(2, 1, 2)
pyplot.imshow(x[ecs].states3d.mean(2).T, extent=x[ecs].extent('xy'), origin="lower", aspect='equal')
fig.tight_layout()
pyplot.show()

With this mod files to describe the simple pump mechanism;

Code: Select all

NEURON {
        SUFFIX pump 
        USEION x READ xi, xo WRITE ix VALENCE 1
        NONSPECIFIC_CURRENT i
        RANGE ix, rate
}

UNITS {
	(mV)	= (millivolt)
	(molar) = (1/liter)
	(mM)	= (millimolar)
	(mA)	= (milliamp)
}


PARAMETER {
    imax    (mA/cm2)
    Kd      (mM)
}

INITIAL {
    imax = 1.0
    Kd = 1.0
}
 
ASSIGNED {
    xo      (mM)
    xi      (mM)
    ix      (mA/cm2)
    i       (mA/cm2)
}
 
BREAKPOINT {
    ix = imax * (xi/(Kd + xi) - xo/(Kd + xo)) 
    i = -ix
}

PS: Beginning in NEURON 7.7, creating a Vector and setting the record can be done in one line, e.g.
t_vec = h.Vector().record(h._ref_t)

But avoid this syntax if your code must run on older versions of NEURON (including the currently released version) as in those versions t_vec would be None and no recording would happen.

[bremen]

Hi.
Thank you.

Currently the easiest way to do this is with mod files describing pump mechanism.

I'll try this approach.

[bremen]

Hi.

After some more testing i have some questions:
1) I have recorded a vector, pointed to the outside name of the substance "self.soma(0.5)._ref_xo", to record the concentration in the ecs. The resulting graph is congruent with my expectation since it records only a part of the entire diffusion.

Is there a way to know or define the location where the recording is taking place?

2) Using this part of you example code..

Code: Select all

pyplot.imshow(x[ecs].states3d.mean(2).T, extent=x[ecs].extent('xy'), origin="lower", aspect='equal')

... and linked to the previous graph, i was expecting to see a "fading yellow core" when the substance concentration started to be lower than the rest of the ecs. At 1000ms the concentration is near zero but the yellow core does not change color.

Maybe i have not clear the exact purpose of states3d. Is there some documentation how to correctly use it?

Best
Stefano

[adamjhn]

1)
The following function will return the location of a segment e.g. 'cell1(0.5)', adapted from rxd.species._xyz;

Code: Select all

def xyz(seg):
    """Return the (x, y, z) coordinate of the center of the segment."""
    sec = seg.sec
    n3d = sec.n3d()
    x3d = [sec.x3d(i) for i in range(n3d)]
    y3d = [sec.y3d(i) for i in range(n3d)]
    z3d = [sec.z3d(i) for i in range(n3d)]
    arc3d = [sec.arc3d(i) for i in range(n3d)]
    return numpy.interp([seg.x * sec.L], arc3d, x3d)[0], numpy.interp([seg.x * sec.L], arc3d, y3d)[0], numpy.interp([seg.x * sec.L], arc3d, z3d)[0]

You can also record from any extracellular location using the node, e.g. for species 'x' on Extracellular space 'ecs' at the origin (0,0,0);

Code: Select all

x[ecs].node_by_location(0,0,0)._ref_concentration

2)
states3d is a 3D numpy array of concentrations in the extracellular grid.
With imshow; “By default, the colormap covers the complete value range of the supplied data.” https://matplotlib.org/api/_as_gen/matp ... mshow.html. So even though the value is lower, it is still the largest value in the data and so is assigned the same color. You can set the range of the colormap with 'vmin' and 'vmax' arguments. You can also show the range with pyplot.colorbar().

[bremen]

Hi.
Thank you very much.

Best regards
Stefano

[bremen]

Hi.

I will soon start working on a particular reaction and, before doing it with a mod file, i tried MultiCompartmentReaction with the lastest code (compiled today).

Code: Select all

ecs_cell2_membrane = rxd.Region([cell2.soma], name='mem', geometry = rxd.ScalableBorder(1, on_cell_surface=False))
pump_ecs_cell2 = rxd.MultiCompartmentReaction(x[ecs] > x[soma_cell2], 10, membrane = ecs_cell2_membrane)

The specie is declared for all the sections involved.

Code: Select all

Traceback (most recent call last):
  File "/usr/local/lib/python2.7/dist-packages/neuron/crxd/rxd.py", line 1239, in _init
    _compile_reactions()
  File "/usr/local/lib/python2.7/dist-packages/neuron/crxd/rxd.py", line 1119, in _compile_reactions
    operator = '+=' if ecs_species_ids_used[species_id][region_id] else '='
IndexError: index 1 is out of bounds for axis 0 with size 1
/nrn77/nrn/x86_64/bin/nrniv: Python Callback failed
 near line 0

Cheers
Stefano