Play Stimulus Vectors Sequentially During Run-time

[ajpc6d]

This code is a proof-of-principle experiment in playing a vector into a point process/distributed mechanism during run-time. (Reason: for extremely complex stimuli, the entire stimulus cannot be contained in memory at once. To circumvent this limitation, the stimulus will be calculated in chunks passed serially to NEURON throughout the simulation.)

Code: Select all

from neuron import h 
import numpy as np
import matplotlib.pyplot as plt
h.load_file("stdrun.hoc")

tvec = h.Vector(np.linspace(0,20,500))
avec = h.Vector(2*np.sin(2*np.pi*tempo))

def set_amp(stim, av, tv):
    stim.dur = 20
    stim.delay = 10
    av.play(stim._ref_amp, tv, True)


st = 100
pw = 3

cell = h.Section(name='cell')
cell.nseg = 5
cell.L = 1e2
cell.diam = 4
cell.insert(h.hh)

stim = h.IClamp(cell(0.5))

v = h.Vector().record(cell(0.5)._ref_v)
t = h.Vector().record(h._ref_t)

h.finitialize(-67)
ss = (set_amp, (stim, avec, tvec))
h.CVode().event(10, ss)
h.continuerun(st)

plt.plot(t,v)
plt.show()

This code does not successfully play anything into the h.IClamp() point process.

In the code below, set_amp()works, so I know a stimulus can be defined at runtime.. but set_amp2() - with the vector play and distributed mechanism - fails.

Code: Select all

# a simple example of the use of h.CVode().event() to pass stimulus data in chunks

from neuron import h 
import numpy as np
import matplotlib.pyplot as plt
h.load_file("stdrun.hoc")

# this function tests setting a source by h.CVode()
# confirmed working
def set_amp_simple(stim,t):
    stim.dur = 20
    stim.delay = 10 + t
    stim.amp = 0.1
   
# this function adds vector play 
# confirmed NOT working
def set_amp2(stim, sa, tv):
    stim.dur = 20
    stim.delay = 10
    sa.play(stim._ref_amp, tv, True)

# define the run time, a time vector, and a stimulus amplitude vector
st = 100
fc=1e3
tvec = h.Vector(np.linspace(0,20,5e3))
stim_amp = h.Vector(2*np.sin(2*np.pi*fc*t))

# define the cell
cell = h.Section(name='cell')
cell.nseg=5
cell.L = 1e2
cell.diam = 4
cell.insert(h.hh)

# define sources
# for set_amp(), each source is turned on once, at a different time, to "layer" stimuli
# for set_amp2(), only one source is needed
stim1 = h.IClamp(cell(0.5))
stim2 = h.IClamp(cell(0.5))
stim_ls = [stim1, stim2]

# vectors to record response
v = h.Vector().record(cell(0.5)._ref_v)
t = h.Vector().record(h._ref_t)
svec = h.Vector()

h.finitialize(-67)
# these lines, for set_amp()
tset = [0,50]
for (T, s) in zip(tset, stim_ls):
    ss = (set_amp_simple, (s, T))
    h.CVode().event(T, ss)
# these lines, for set_amp2()
# ss = (set_amp2, (stim1, stim_amp, tvec))
# h.CVode().event(10, ss)
h.continuerun(st)

plt.plot(t,v)
plt.show()

[ted]

for extremely complex stimuli, the entire stimulus cannot be contained in memory at once

This is generally a sign that the stimulus should be computed on the fly, especially when it can be implemented nicely with NMODL. Such is the case with many mathematical formulas such as

2*np.sin(2*np.pi*tempo)

which is mathematically equivalent to the example provided by Izap, which we have discussed elsewhere.

If you have some other function in mind that you would prefer to discuss privately, feel free to email me.

[ajpc6d]

Done - thanks so much