To run the code below:
SHIFT+ENTER
on your keyboard or press the play button
() in the toolbar above.Feel free to create new cells using the plus button
(), or pressing SHIFT+ENTER
while this cell
is selected.
This is an interactive verison of the idealized model of the smooth pursuit reflex. This version does not explain the model itself, but shows how Brian's "runtime mode" can be used to interact with a running simulation. In this mode, the generated code based on the model descriptions is seemlessly integrated with the Python environment and can execute arbitrary Python code at any point during the simulation via a specially annotated function, called a "network operation".
For a non-interactive version of this example which generates the article's figure see this notebook.
This version uses the plotly library for consistency with the other examples. However, updating an existing figure is slower than with the matplotlib library, so it is not an ideal solution for this example. For a version based on matplotlib, see this notebook
# Needs ipywidgets
import ipywidgets as widgets
import threading
from brian2 import *
The model itself (mostly identical to the non-interactive example, except that some of the constants are included as parameters in the equation and can therefore change during the simulation):
alpha = (1/(50*ms))**2 # characteristic relaxation time is 50 ms
beta = 1/(50*ms) # friction parameter
eqs_eye = '''
dx/dt = velocity : 1
dvelocity/dt = alpha*(x0-x)-beta*velocity : 1/second
dx0/dt = -x0/tau_muscle : 1
dx_object/dt = (noise - x_object)/tau_object: 1
dnoise/dt = -noise/tau_object + tau_object**-0.5*xi : 1
tau_object : second
tau_muscle : second
'''
eye = NeuronGroup(1, model=eqs_eye, method='euler')
taum = 20*ms
motoneurons = NeuronGroup(2, model= 'dv/dt = -v/taum : 1', threshold = 'v>1',
reset = 'v=0', refractory = 5*ms, method='exact')
motosynapses = Synapses(motoneurons, eye, model = 'w : 1', on_pre = 'x0+=w')
motosynapses.connect() # connects all motoneurons to the eye
motosynapses.w = [-0.5,0.5]
N = 20
width = 2./N # width of receptive field
gain = 4.
eqs_retina = '''
I = gain*exp(-((x_object-x_eye-x_neuron)/width)**2) : 1
x_neuron : 1 (constant)
x_object : 1 (linked) # position of the object
x_eye : 1 (linked) # position of the eye
dv/dt = (I-(1+gs)*v)/taum : 1
gs : 1 # total synaptic conductance
'''
retina = NeuronGroup(N, model = eqs_retina, threshold = 'v>1', reset = 'v=0', method='exact')
retina.v = 'rand()'
retina.x_eye = linked_var(eye, 'x')
retina.x_object = linked_var(eye, 'x_object')
retina.x_neuron = '-1.0 + 2.0*i/(N-1)'
sensorimotor_synapses = Synapses(retina, motoneurons, model = 'w : 1 (constant)', on_pre = 'v+=w')
sensorimotor_synapses.connect(j = 'int(x_neuron_pre > 0)')
sensorimotor_synapses.w = '20*abs(x_neuron_pre)/N_pre'
M = StateMonitor(eye, ('x', 'x0', 'x_object'), record = True)
S_retina = SpikeMonitor(retina)
S_motoneurons = SpikeMonitor(motoneurons)
We create an empty plot that will be updated during the run:
from plotly import tools
from plotly.offline import iplot, init_notebook_mode
import plotly.graph_objs as go
init_notebook_mode(connected=True)
fig = tools.make_subplots(3, 1, specs=[[{'rowspan': 2}], [None], [{}]],
shared_xaxes=True, print_grid=False)
retina_spikes = go.Scattergl(x=[],
y=[],
marker={'symbol': 'line-ns', 'line': {'width': 1, 'color':'black'}},
mode='markers',
name='retina',
showlegend=False)
fig.append_trace(retina_spikes, 1, 1)
left_spikes = go.Scattergl(x=[],
y=[],
marker={'symbol': 'line-ns', 'line': {'width': 1, 'color':'#1f77b4'},
'color':'#1f77b4'},
mode='markers',
name='left motoneuron',
showlegend=False)
fig.append_trace(left_spikes, 1, 1)
right_spikes = go.Scattergl(x=[],
y=[],
marker={'symbol': 'line-ns', 'line': {'width': 1, 'color':'#ff7f03'},
'color':'#ff7f03'},
mode='markers',
name='right motoneuron',
showlegend=False)
fig.append_trace(right_spikes, 1, 1)
eye_position = go.Scattergl(x=[],
y=[],
mode='lines',
line={'color': 'black'},
name='eye')
fig.append_trace(eye_position, 3, 1)
object_position = go.Scattergl(x=[],
y=[],
mode='lines',
line={'color': '#2ca02c'},
name='object')
fig.append_trace(object_position, 3, 1)
fig['layout'].update(xaxis1={'showline': False,
'zeroline': False,
'title': 'time (in s)',
'range': (0, 10),
'ticktext': ['-10', '-8', '-6', '-4', '-2', '0'],
'tickvals': [0, 2, 4, 6, 8, 10]},
yaxis1={'title': 'neuron index',
'showticklabels': False,
'showline': True,
'range': (0, N+2)},
yaxis2={'tickmode': 'array',
'ticktext': ['left', 'right'],
'tickvals': [-1, 1],
'range': [-1.05, 1.05],
'zeroline': True,
'showline': True}
)
fig_widget = go.FigureWidget(fig)
retina_spikes, left_spikes, right_spikes, eye_position, object_position = fig_widget.data
We now create interactive widgets that the user can use to start/stop the simulation, as well as for setting certain simulation parameters.
time_label = widgets.Label(value='Time: 0 s')
start_stop_button = widgets.Button(tooltip='Start simulation', icon='play')
tau_obj_slider = widgets.FloatSlider(orientation='horizontal', description='tau_object',
value=500, min=100, max=1000)
tau_muscle_slider = widgets.FloatSlider(orientation='horizontal', description='tau_muscle',
value=20, min=5, max=100)
weight_slider = widgets.FloatSlider(orientation='horizontal', description='w_muscle',
value=0.5, min=0, max=2)
sliders = widgets.VBox([widgets.HBox([time_label, start_stop_button]),
tau_obj_slider, tau_muscle_slider, weight_slider])
layout = widgets.HBox([fig_widget, sliders])
We interact with the running simulation via a "network operation", a Python function that will be regularly called by Brian during the simulation run (here, every 100ms of biological time). This function can access arbitrary attributes of the model to get or set their values. We use this here to 1) update the plot with the data from the last second and 2) set parameters of the model to the values requested by the user.
should_stop = False
@network_operation(dt=500*ms)
def plot_output(t):
with fig_widget.batch_update():
cutoff = (t - 10*second)
# Plot the data of the last 10 seconds
indices = S_retina.t > cutoff
retina_spikes.x = (S_retina.t[indices] - cutoff)/second
retina_spikes.y = S_retina.i[indices]
motoneuron_trains = S_motoneurons.spike_trains()
to_plot = motoneuron_trains[0][motoneuron_trains[0] > cutoff]
left_spikes.x = (to_plot - cutoff)/second
left_spikes.y = np.ones(len(to_plot))*N
to_plot = motoneuron_trains[1][motoneuron_trains[1] > cutoff]
right_spikes.x = (to_plot - cutoff)/second
right_spikes.y = np.ones(len(to_plot))*(N+1)
indices = M.t > cutoff
eye_position.x = (M.t[indices] - cutoff)/second
eye_position.y = M.x[0][indices]
object_position.x = (M.t[indices] - cutoff)/second
object_position.y = M.x_object[0][indices]
time_label.value = 'Time: {:.1f}s'.format(float(t[:]))
# Set the simulation parameters according to user settings
eye.tau_object = tau_obj_slider.value*ms
eye.tau_muscle = tau_muscle_slider.value*ms
motosynapses.w = [-weight_slider.value, weight_slider.value]
if should_stop:
net.stop()
We store the model and the "network operation" in a Network
object, and store its current state to allow for repeated execution.
net = Network(collect())
net.store()
We now define two helper functions used to start/stop simulations. The actual simulation will be run in a background thread so that the user interface stays reactive while the simulation is running:
def do_run(runtime):
net.restore()
net.run(runtime)
running = False
def button_pressed(b):
global running
global should_stop
if running:
should_stop = True
running = False
start_stop_button.tooltip = 'Start simulation'
start_stop_button.icon = 'play'
else:
should_stop = False
running = True
time_label.value = 'starting...'
start_stop_button.tooltip = 'Stop simulation'
start_stop_button.icon = 'stop'
thread = threading.Thread(target=do_run, args=(100*second, ))
thread.start()
start_stop_button.on_click(button_pressed)
We are now ready to display the plot and user interface, which can then be used to start the simulation and interact with the simulation parameters:
display(layout)