Store/Restore Stochastic Disruption on Start up

[cwalsh27]

Description of problem

I’m using the store/restore system in my model, but it’s not performing as I expected. I run the model for a 50-second spin-up period, and then store the network. I then restore the network and run it for a 300 ms experimental period. However, when I first restore the network, the PoissonInput objects that I’m using for noise in the network seem to spike before settling down again. Is this expected behavior, or is there something that I’m doing incorrectly that is causing the jump?

I’ll also add that I am performing these tests in a Jupyter notebook, and the store/restore are in different cells, but I don’t imagine that to be the cause?

Minimal code to reproduce problem

# network definition 

Thal.I_app = 0*pA
L4e.I_app = 0*pA
L4i.I_app = 0*pA
Ret.I_app = 0*pA
Spike_t = SpikeMonitor(Thal)
Spike_e = SpikeMonitor(L4e)
Spike_i = SpikeMonitor(L4i)
Spike_r = SpikeMonitor(Ret)

run(Init_time, report='stdout')
store()

# a few steps dedicated to establishing directories for output data

restore()

del Spike_t, Spike_e, Spike_i, Spike_r
Spike_t = SpikeMonitor(Thal)
Spike_e = SpikeMonitor(L4e)
Spike_i = SpikeMonitor(L4i)
Spike_r = SpikeMonitor(Ret)
setParam('Istim', str(i) + '*pA', params)
setParam('StimCx', str(s), params)

Thal.I_app = Istim * (1-StimCx)
L4e.I_app = Istim * StimCx
L4i.I_app = Istim * StimCx

run(run_time)

What you have aready tried

I’ve implemented a buffer window, allowing the model to run an extra 1.5 s before I reactivate the spike counters. That way, even though the spike still happens, I get a clean reading from a stable network. I’m more curious about whether or not I’m doing something wrong, or if this is the intended functionality of the store()/restore() commands.

Actual output

Below is a plot of the network activity at the end of the spin-up period and the beginning of the run-time period. The spin-up is on the left and the run-time is on the right, so the two plots can be viewed as one time axis side by side. As shown, at the beginning of the run-time period, a level of disruption occurs when the network first reactivates.

[mstimberg]

Hi @cwalsh27 Thanks for the report! Unfortunately, your example is missing the network definition part, so I cannot run it. I tried a simple example with PoissonInput (and PoissonGroup), but did not see anything off ­– I guess you are doing something different, though. Here’s my code for reference:

from brian2 import *

group = NeuronGroup(1, """
                       dx/dt = -x / (10*ms) : 1
                       dy/dt = -y / (10*ms) : 1  
                       """)
poisson_1 = PoissonGroup(100, rates=100*Hz)
syn_poisson = Synapses(poisson_1, group, on_pre="x += 1")
syn_poisson.connect()
poisson_2 = PoissonInput(group, "y", 100, 100*Hz, weight=1)

mon = StateMonitor(group, ["x", "y"], record=0)
run(1*second)
store()
restore()
run(500*ms)
plt.plot(mon.t/ms, mon.x[0])
plt.plot(mon.t/ms, mon.y[0])
plt.axvline(1000, color="darkgray", zorder=-1)
plt.show()

which gives