example_2_eye_movements_interactive_matplotlib.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Quickstart\n",
    "To run the code below:\n",
    "\n",
    "1. Click on the cell to select it.\n",
    "2. Press `SHIFT+ENTER` on your keyboard or press the play button\n",
    "   (<button class='fa fa-play icon-play btn btn-xs btn-default'></button>) in the toolbar above.\n",
    "\n",
    "Feel free to create new cells using the plus button\n",
    "(<button class='fa fa-plus icon-plus btn btn-xs btn-default'></button>), or pressing `SHIFT+ENTER` while this cell\n",
    "is selected."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Example 2 (Smooth pursuit eye movements) – interactive version based on *matplotlib*\n",
    "\n",
    "\n",
    "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\".\n",
    "\n",
    "For a non-interactive version of this example which generates the article's figure see [this notebook](example_2_eye_movements.ipynb).\n",
    "\n",
    "This notebook is based on *matplotlib* and *ipympl*, which enables quick updates of the plot in real-time. For a version based on *plotly* (as the other, non-interactive examples), see [this notebook](example_2_eye_movements_interactive.ipynb)."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Needs ipywidgets and ipympl\n",
    "%matplotlib widget\n",
    "import ipywidgets as widgets\n",
    "import threading\n",
    "from brian2 import *\n",
    "\n",
    "plt.ioff()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The model itself (mostly identical to the [non-interactive example](example_2_eye_movements.ipynb), except that some of the constants are included as parameters in the equation and can therefore change during the simulation):"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "alpha = (1/(50*ms))**2 # characteristic relaxation time is 50 ms\n",
    "beta = 1/(50*ms) # friction parameter\n",
    "\n",
    "eqs_eye = '''\n",
    "dx/dt = velocity : 1\n",
    "dvelocity/dt = alpha*(x0-x)-beta*velocity : 1/second\n",
    "dx0/dt = -x0/tau_muscle : 1\n",
    "dx_object/dt = (noise - x_object)/tau_object:  1\n",
    "dnoise/dt = -noise/tau_object + tau_object**-0.5*xi : 1\n",
    "tau_object : second\n",
    "tau_muscle : second\n",
    "'''\n",
    "eye = NeuronGroup(1, model=eqs_eye, method='euler')\n",
    "taum = 20*ms\n",
    "motoneurons = NeuronGroup(2, model= 'dv/dt = -v/taum : 1', threshold = 'v>1',\n",
    "                          reset = 'v=0', refractory = 5*ms, method='exact')\n",
    "motosynapses = Synapses(motoneurons, eye, model = 'w : 1', on_pre = 'x0+=w')\n",
    "motosynapses.connect() # connects all motoneurons to the eye\n",
    "motosynapses.w = [-0.5,0.5]\n",
    "N = 20\n",
    "width = 2./N # width of receptive field\n",
    "gain = 4.\n",
    "eqs_retina = '''\n",
    "I = gain*exp(-((x_object-x_eye-x_neuron)/width)**2) : 1\n",
    "x_neuron : 1 (constant)\n",
    "x_object : 1 (linked) # position of the object\n",
    "x_eye : 1 (linked) # position of the eye\n",
    "dv/dt = (I-(1+gs)*v)/taum : 1\n",
    "gs : 1 # total synaptic conductance\n",
    "'''\n",
    "retina = NeuronGroup(N, model = eqs_retina, threshold = 'v>1', reset = 'v=0', method='exact')\n",
    "retina.v = 'rand()'\n",
    "retina.x_eye = linked_var(eye, 'x')\n",
    "retina.x_object = linked_var(eye, 'x_object')\n",
    "retina.x_neuron = '-1.0 + 2.0*i/(N-1)'\n",
    "sensorimotor_synapses = Synapses(retina, motoneurons, model = 'w : 1 (constant)', on_pre = 'v+=w')\n",
    "sensorimotor_synapses.connect(j = 'int(x_neuron_pre > 0)')\n",
    "sensorimotor_synapses.w = '20*abs(x_neuron_pre)/N_pre'\n",
    "M = StateMonitor(eye, ('x', 'x0', 'x_object'), record = True)\n",
    "S_retina = SpikeMonitor(retina)\n",
    "S_motoneurons = SpikeMonitor(motoneurons)\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We create an empty plot that will be updated during the run:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Plot preparation\n",
    "fig, (ax_spikes, ax_position) = plt.subplots(2, 1, gridspec_kw={'height_ratios': (2, 1)}, sharex=True)\n",
    "h_retina = ax_spikes.plot([], [], '|k', markeredgecolor='k', label='retina')[0]\n",
    "h_left = ax_spikes.plot([], [], '|', color='C0', markeredgecolor='C0', label='left motoneuron')[0]\n",
    "h_right = ax_spikes.plot([], [], '|', color='C1', markeredgecolor='C1', label='right motoneuron')[0]\n",
    "ax_spikes.set(yticks=[], ylabel='neuron index', xticks=[], xlim=(0, 10), ylim=(0, 22))\n",
    "ax_spikes.spines['bottom'].set_visible(False)\n",
    "\n",
    "ax_position.axhline(0, color='gray')\n",
    "h_eye = ax_position.plot([], [], 'k', label='eye')[0]\n",
    "h_object = ax_position.plot([], [], color='C2', label='object')[0]\n",
    "ax_position.set(yticks=[-1, 1], yticklabels=['left', 'right'], xlabel='time (s)',\n",
    "                xticks=np.arange(11, 2), xticklabels=np.arange(11, 2)-10,\n",
    "                xlim=(0, 10), ylim=(-1, 1))\n",
    "ax_position.legend(loc='upper right', bbox_to_anchor=(1.0, 2.0));"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We now create interactive widgets that the user can use to start/stop the simulation, as well as for setting certain simulation parameters."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "time_label = widgets.Label(value='Time: 0 s')\n",
    "start_stop_button = widgets.Button(tooltip='Start simulation', icon='play')\n",
    "\n",
    "tau_obj_slider = widgets.FloatSlider(orientation='horizontal', description='tau_object',\n",
    "                                     value=500, min=100, max=1000)\n",
    "tau_muscle_slider = widgets.FloatSlider(orientation='horizontal', description='tau_muscle',\n",
    "                                        value=20, min=5, max=100)\n",
    "weight_slider = widgets.FloatSlider(orientation='horizontal', description='w_muscle',\n",
    "                                    value=0.5, min=0, max=2)\n",
    "sliders = widgets.VBox([widgets.HBox([time_label, start_stop_button]),\n",
    "                        tau_obj_slider, tau_muscle_slider, weight_slider])\n",
    "layout = widgets.HBox([fig.canvas, sliders])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "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."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "should_stop = False\n",
    "@network_operation(dt=100*ms)\n",
    "def plot_output(t):\n",
    "    cutoff = (t - 10*second)\n",
    "    # Plot the data of the last 10 seconds\n",
    "    indices = S_retina.t > cutoff\n",
    "    h_retina.set_data((S_retina.t[indices] - cutoff)/second, S_retina.i[indices])\n",
    "    motoneuron_trains = S_motoneurons.spike_trains()\n",
    "    to_plot = motoneuron_trains[0][motoneuron_trains[0] > cutoff]\n",
    "    h_left.set_data((to_plot - cutoff)/second, np.ones(len(to_plot))*N)\n",
    "    to_plot = motoneuron_trains[1][motoneuron_trains[1] > cutoff]\n",
    "    h_right.set_data((to_plot - cutoff)/second, np.ones(len(to_plot))*(N+1))\n",
    "    indices = M.t > cutoff\n",
    "    h_eye.set_data((M.t[indices] - cutoff)/second, M.x[0][indices])\n",
    "    h_object.set_data((M.t[indices] - cutoff)/second, M.x_object[0][indices])\n",
    "    fig.canvas.draw_idle()\n",
    "    time_label.value = 'Time: {:.1f}s'.format(float(t[:]))\n",
    "    # Set the simulation parameters according to user settings\n",
    "    eye.tau_object = tau_obj_slider.value*ms\n",
    "    eye.tau_muscle = tau_muscle_slider.value*ms\n",
    "    motosynapses.w = [-weight_slider.value, weight_slider.value]\n",
    "    if should_stop:\n",
    "        net.stop()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We store the model and the \"network operation\" in a `Network` object, and store its current state to allow for repeated execution."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "net = Network(collect())\n",
    "net.store()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "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:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def do_run(runtime):\n",
    "    net.restore()\n",
    "    net.run(runtime)\n",
    "running = False\n",
    "def button_pressed(b):\n",
    "    global running\n",
    "    global should_stop\n",
    "    if running:\n",
    "        should_stop = True\n",
    "        running = False\n",
    "        start_stop_button.tooltip = 'Start simulation'\n",
    "        start_stop_button.icon = 'play'\n",
    "    else:\n",
    "        should_stop = False\n",
    "        running = True\n",
    "        time_label.value = 'starting...'\n",
    "        start_stop_button.tooltip = 'Stop simulation'\n",
    "        start_stop_button.icon = 'stop'\n",
    "        thread = threading.Thread(target=do_run, args=(100*second, ))\n",
    "        thread.start()\n",
    "start_stop_button.on_click(button_pressed)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "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:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "display(layout)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.6.7"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}