Turning current on and off and connecting it into a subpopulation of neurons

Hi there,

I am trying to define a code which turns an input current (Ornstein-Uhlenbeck process) on and off at specific times, and connect it to a small specific subpopulation of neuron group E with a specific weight:

E=NeuronGroup(…)
I=NeuronGroup(…)

myclock=Clock(dt=1*ms)
@network_operation(myclock)
def update_input():
if myclock.t >= stim_on and myclock.t < stim_off:
I1 = i1[:,int( (myclock.t - stim_on) / (1 * ms))] * amp
I2 = i2[:,int( (myclock.t - stim_on) / (1 * ms))] * amp
else:
I1 = 0 * amp
I2 = 0 * amp

net = Network(E, I)
net.run(runtime)

Now, if I am not mistaken, I1 and I2 are formed and I am going to insert these currents to two distinct subpopulations of E neuron groups, lets say E1 and E2.
How should I do that? Defining different subpopulation via different neuron groups? or …
To the best of my knowledge, in Brian1 that was simply done by E1.I =I1,etc.

Regards,
Hedyeh

Hi @Hedyeh,
you can do quite a bit more in equations than what Brian 1 supported, so you don’t have to use a (relatively inefficient) network_operation here.
To switch some input only on during a certain period, you have two options:

  1. Use several run statements. E.g. in your equations, use something like int(input_active)*I, i.e. multiply the input by 0 or 1, then set the flag accordingly between runs:
input_active = False
run(stim_on)
input_active = True
run(stim_off - stim_on)
input_active = False
run(runtime - stim_off)
  1. multiply the input current in the equations with an expression of time. I.e., instead of using a separate input_active variable as in the example above, directly use int(t >= stim_on and t<stim_off)*I.

Regarding different inputs to different cells, a TimedArray can be 2D, and you then use it as an argument of time and neuron index:

I = TimedArray(..., dt=1*ms)
eqs = '''
...
... = ... * I(t, i) ...
'''

In the example above, every neuron gets its own input (which could also be just zeros for some of them), but this can be much more flexible, if needed. For example, you could use I(t, i//(N/2)) so that the first half of the neurons gets one type of input and the second half gets another type. For full flexibility, you could also define a number of input types (i.e. columns of the matrix provided to TimedArray), and then use a specific variable that denotes the type for each neuron:

eqs = '''
input_type : integer(constant)
... = ...*I(t, input_type) ...
'''
E = NeuronGroup(1000, eqs, ...)
# input_type is 0 by default
E1 = E[100:200]
E1.input_type = 1
E2 = E[200:300] 
E2.input_type = 2

I hope that makes the idea clear. Best
Marcel

Dear Marcel,
Thanks a lot for the detailed reply, I just produced the timed array as:
I = TimedArray(numpy.r_[numpy.zeros(200), numpy.ones(200), numpy.zeros(200)], dt = 1 * ms)
Sorry, I could not figure out how it can be extended to a two dimensional one to include neuron indices as well? and how should it be included in neuron equations? And finally how this current can carry two different values for two different subopulations? I need this current to be injected to E1 subpop with value i1 for the middle 200 ms, and to E2 subpop with value i2 for the same time interval.
E1 = E[N_non:N_non+N_sub]
E2 = E[N_non+N_sub: N_non+2*N_sub]

And i1 and i2 are filtered Ornstein-Uhlenbeck process with the length of 200 ms.

eqs_E = ‘’’
dv/dt = (I - v) / C : volt
I : amp
‘’’
E = NeuronGgroup(N_E, eqs_E, …)
E.I = 0

And would also be so kind as to let me know more about network_operation in this case? What needs to be included in net = network() in case I need to resort to network_operation here? All neurons and synapses?

Thank you so much in advance.

Regards,
Hedyeh

Hi @Hedyeh.
The general usage of TimedArray is explained in the docs: Input stimuli — Brian 2 2.5.1 documentation
I think the clearest/simple approach for your issue would be to have an array with three columns, one for the population that should not get any inputs (so zero everywhere), and one each for the populations. So something like this (for simplicity I am using only 1/100th of the values here):

i1 = numpy.random.randn(2)
i2 = numpy.random.randn(2)
values = numpy.vstack([numpy.zeros(6),
                       numpy.r_[numpy.zeros(2), i1, numpy.zeros(2)],
                       numpy.r_[numpy.zeros(2), i2, numpy.zeros(2)]]).T

In this example, values now look something like this:

array([[ 0.        ,  0.        ,  0.        ],
       [ 0.        ,  0.        ,  0.        ],
       [ 0.        , -1.3455429 , -0.88642223],
       [ 0.        ,  1.99887995,  0.0135776 ],
       [ 0.        ,  0.        ,  0.        ],
       [ 0.        ,  0.        ,  0.        ]])

You then use the approach that I explained earlier, where you have a special constant that defines for each neuron what input it should take. Note that you should not define I as a variable of your model – a TimedArray in brian2 is used like a function:

I = TimedArray(values*nA, dt=1*ms)
eqs_E = '''
dv/dt = I(t, input_type)/ C : volt
input_type : integer (constant)
'''
E = NeuronGroup(N_E, eqs_E, …)
E1 = E[...]
E2 = E[...]
E1.input_type = 1
E2.input_type = 2
# input_type is 0, i.e. "first column", for all other neurons

Dear Marcel,
Many thanks for your help. I was able to produce different currents for each neuron across each subpopulation (I mean not the same current for all neuron within the subpopulation). Then, I need to have an array for values with dimensions (time*N_E1+N_E2). I just copy the piece of code here in case of help to others.
Best,
Hedyeh

i1 and i2 random numbers with size (stimulus-duration*subpopulation-size) to be injected to neurons in subpop1 and subpop2

values = np.hstack([np.r_[np.zeros((int(s_on), N_E1)), i1,
np.zeros((int(r_ti-s_of), N_E1))],
np.r_[np.zeros((int(s_on), N_E2)), i2,
np.zeros((int(r_ti-s_of), N_E2))]])

stimulus = TimedArray(values * nA, dt=1*ms)

Different input for each subgroup

pop_E1.input_type = np.arange(0,N_E1)
pop_E2.input_type = np.arange(N_E1,2*N_E1)