Brian 2

  1. What it is

  2. Installation

  1. Brian2 Basics

  1. Spiking Neuron Models in Brian2

  2. Literature Examples

  3. Strengths and Weaknesses

  4. Resources
  5. references




What it is

  • Python package for simulating biologically detailed networks of spiking neurons

  • Can simulate neuronal activity at multiple scales:

    • Dynamics between multiple interconnected populations
    • Population-level activity
    • Detailed activity of small networks
    • Spatially accurate biophysical characteristics of individual neurons
  • Flexible - the level of detail used can be adapted to the needs of the experiment and the researcher.


Installation

Brian2 is simple to install

With conda:
 conda create -n <your_env>
 activate <your_env>
 conda install -c conda-forge brian2
 conda install matplotlib nose ipython jupyter

The brian2tools package provides some additional functions for graphing collected data and neuron morphologies

conda install -c brian_team brian2tools

Importing brian2

brian2 works best when the whole module is imported into the namespace. This import includes matplotlib and numpy functions configured to work with brian2's unit system.

from brian2 import *

to prefix numpy and matplotlib functions while keeping their brian2 versions intact, use the import structure below

from brian2 import *
import brian2.numpy_ as np
import matplotlib.pyplot as plt
Note :

the full brian2 import brings a function named input that interferes with the built-in python function of the same name. The function below can be used to accept user input without interfering with brian2

def input_wrap(inpt):
    '''Wraps built-in input() function to fix name conflict with brian2's
    input() function. inpt: type(str)'''
    user_in = __builtins__.input(inpt)
    return u_in




Brian2 Basics

Units

Brian2 uses the SI unit system, and calculations return dimensionally consistent results. For example in Ohm's law $V = IR$ , $I$ has units of amps and $R$ has units of Ohms

With brian2's unit system, this would be written as:

I = 2*amp
R = 550*ohm
V = I * R

When evaluated, the voltage is returned with the correct unit and scale

1.1 kV

Key components of brian2 code

The following elements form the base of most code written for brian2 simulations.

Equations :

eqs = Equations('''
dv/dt = (I - v)/tau : 1 <any flags>
I : <unit I>
tau : <unit tau>
''')
  • define the activity of models with differential equations, and the state variables within models
  • pass to the model= parameter of NeuronGroup and Synapses objects



NeuronGroup :

G = NeuronGroup(N=10,model=eqs,threshold='v>0.9',reset='v=0',refractory=3*ms,method='euler')
  • creates a population of N neurons that use the same model equations
  • define spiking behavior for that population; method='euler' defines the type of numerical integration to use
  • Note - object variables can be accessed without their associated units by adding an underscore to the end, e.g. "G.v_"


Synapses :

S = Synapses(G, G1, model=''' 
w : 1
p : 1''', on_pre='v_post += w * (rand()<p)', on_post=None)

S.w = 0.4
S.p = 0.8
  • define dynamics of synaptic variables and connections between neurons in the presynaptic (G) and postsynaptic (G1) NeuronGroups
    • here, synaptic weight S.w and the chance it will be applied S.p are respectively set at 0.4 and 0.8 for the whole synapse
  • on_pre and on_post define what happens when pre and postsynaptic spikes fire
  • synaptic connections can be instantiated by index
    • S.connect(i=1, j=4) connects presynaptic neuron 1 to postsynaptic neuron 4
  • or probability
    • S.connect(p=0.5) has a 50% chance of connecting any two neurons
    • S.connect(p='abs(i-j) /(1.0*n_neurons)') yields probability from an equation
  • or condition
    • S.connect(condition='abs(i-j)>2 and abs(i-j) < 5') only connects neurons within certain indices of eachother


StateMonitor :

state_mon = StateMonitor(source=G, variables='v', record=range(10)
  • records values of the state variables in variables from the indices of object source specified in record for NeuronGroups and Synapses


SpikeMonitor :

spike_mon = SpikeMonitor(source=G, record=True)
  • records the spike times from the indices of object source (typically a NeuronGroup) specified in record
  • here, record=True means that spike times will be recorded from all indices in G


Commands : start_scope() and run()

start_scope()
  • stops objects preceding the call from being included in the next run 
run(300*ms)
  • runs the simulation for the specified length of time

Network Operations

@network_operation(dt=10*ms)

def reset_v_every_dt():
    '''set each neuron in G's voltage to 0 after each timstep'''
    for idx in range(n_neurons):
        G.v[idx] = 0
  • allows user-written python code to be executed at every dt in the run
  • great for defining functions that wouldn't otherwise fit well in brian2




Examples

Input-Frequency curves with noise

  • Input-Frequency curves describe the relationship between the input stimulus a neuron receives and its firing rate
  • biological systems are noisy - brian2 has a built-in variable xi (drawn from $X \sim N(0,1)$) to introduce noise into a system
  • The neurons in each group below receive an input current that scales with their index

Spiking Neuron Models in Brian2


Neurobiology of Common Model Parameters

  • Neurons spike when their membrane voltage ($v$) changes enough to reach a certain threshold

  • This change in voltage represents the build-up of opposite electrical charges on the interior and exterior of the membrane; these charges are carried by ionic currents ($I$) flowing through open ion channels in the membrane.

  • Current flow across a membrane depends on how many ion channels are present in that membrane, and how many of those channels are open for current to pass through. The obstruction to current flow caused by these characteristics is quantified in the membrane Resistance $R$ (or $R_M$)

  • When current stops flowing across a membrane, built-up charges will dissipate over time and the membrane voltage will return to $v_{rest}$. The rate at which the membrane charges and discharges is captured in time constant tau ($\tau$)

Models in Brian2 are defined by differential equations

Brian2 can implement models that span a wide range of biological complexity. Though increased realism is preferable in many experiments, more realistic models take more processing power to run and can become impractical for large simulations. Leaky Integrate and Fire (LIF) is a model that trades some biological for the ability to run at much larger scales whereas Hodgkin - Huxley (HH) models make the opposite trade.

Leaky Integrate and Fire

$$\frac{dv}{dt} = \frac{v_{rest} -v(t) + RI(t)}{\tau}$$


  • simple
  • concerned only with spiking - good for simulating large populations
  • neuron spikes when membrane voltage ($v$) reaches a predefined threshold, and $v$ resets
  • computationally light

  • can be simplified further: $$\frac{dv}{dt} = \frac{1-v(t)}{\tau}$$



    Hodgkin - Huxley Model

    $$\frac{dv}{dt} = \frac{I_{L} + I_{K} + I_{Na}}{C_{m}}$$

  • complex
  • accounts for factors like neuron size, ionic conductance, driving force, channel state, etc.

  • computationally heavy - $I_{L}$, $I_{Na}$, $I_{K}$, and the variables $m$, $n$, and $h$ each have their own equation

    • Leak Current $I_L$: $$ I_{L} = g_{L} (E_{L} - v)$$

    • Sodium Current $I_{Na}$: $$ I_{Na} = g_{na} m^{3} h (v - E_{Na}) $$

    • Potassium Current $I_K$: $$ I_{K} = g_{K} n^{4} (v - E_{K}) $$

    Expanded HH equation:
$$\frac{dv}{dt} = g_{L} (E_{L} - v) - g_{na} m^{3} h (v - E_{Na}) - g_{K} n^{4} (v - E_{K}) + I) C_{m}^{-1} $$

Differential equations are passed in as strings

The equations that define model activity are passed as multi-line strings to Equations objects, and each line is separated into two sides by a colon :. The left side defines how the parameters in that line evolve over time, and the right side defines their units and any conditions for the parameters to follow

Leaky Integrate and Fire

eqs_lif_simple = Equations('''
dv/dt = (1-v)/tau : 1 (unless refractory)
''')

eqs_lif = Equations('''
dv/dt = (-v + R * I)/tau : volt (unless refractory)
I : amp
tau : second
''')
  • one differential equation defines the model

  • in the simplification, : 1 indicates that the $v$ variable is unitless and $\tau$ denotes an external variable

  • in the complex version, volt indicates that $v$ is now defined with units, I : amp defines an state variable I with the unit amps, andtau : second defines a state variable $\tau$ with the unit s

    • (unless refractory) specifies a condition in which the model will not follow the equation while in an externally defined refractory state.

Hodgkin - Huxley

eqs_HH = Equations('''
dv/dt = (gl*(El-v) - g_na*(m*m*m)*h*(v-ENa) - g_kd*(n*n*n*n)*(v-EK) + I)/Cm : volt
dm/dt = 0.32*(mV**-1)*(13.*mV-v+VT)/
(exp((13.*mV-v+VT)/(4.*mV))-1.)/ms*(1-m)-0.28*(mV**-1)*(v-VT-40.*mV)/
(exp((v-VT-40.*mV)/(5.*mV))-1.)/ms*m : 1
dn/dt = 0.032*(mV**-1)*(15.*mV-v+VT)/
(exp((15.*mV-v+VT)/(5.*mV))-1.)/ms*(1.-n)-.5*exp((10.*mV-v+VT)/(40.*mV))/ms*n : 1
dh/dt = 0.128*exp((17.*mV-v+VT)/(18.*mV))/ms*(1.-h)-4./(1+exp((40.*mV-v+VT)/(5.*mV)))/ms*h : 1
I : amp
tau : second
''')
  • The equations required to define the Hodgkin - Huxley model are more involved and require additional variables like El to be defined outside the string
  • Equations for the evolution of $m$, $n$, $h$, are defined on their own individual lines

Example: Integrate and Fire neurons driven by semi- random step protocol

The parameters in the examples below can be adjusted to see their effects on spiking activity

GridSpec(6, 4) gs

Example: Hodgkin-Huxley Neurons driven by semi-random step protocol

As above, the parameters in the examples below can be adjusted - Hodgkin - Huxley models can model conditions in much more detail thanks to their level of realism. Toxins and disease states can affect very specific parameters - Tetrodotoxin from pufferfish, for example, would reduce the value of gNa

GridSpec(5, 4) gs

Spatial Representation in Brian2

  • brian2 can simulate spatially accurate neuronal activity with comparable detail to NEURON. Bono and Clopath (2017) used brian2 to evaluate mechanisms of synaptic plasticity at different locations along the basal dendrites of detailed neuron models
    • A comparison between NEURON and brian2 is available in Resources
    • Bono & Clopath (2017) is available in References.
Out[13]:
(-235.32676538123317, 308.9001608273013, -234.6478762380449, 514.9256104152932)

  • Left : The shaded red and blue boxes represent proximal and distal locations on basal dendrites, respectively

    Right : the brian2tools function plot_morphology allows the neuron's structure to be visualized

Out[14]:
(-0.5, 1370.5, 758.5, -0.5)
  • code for this model can be found on ModelDB

Example: Modelling eye movement during smooth pursuit - Stimberg et al., 2019


  • Because models in brian2 are based off of user-defined equations and variables, they can implement diverse phenomena.
  • Stimberg et al. show how brian2 can be used to model eye movements while following a target. Synapses in the code below are used to model the relationships between these groups of neurons
To create this figure, NeuronGroups were used to model:
  • the position and movement of the eye itself : eqs_eye and the eye NeuronGroup
  • motor neurons that control the activity of opposing extra-ocular muscles responsible for moving the eye : motoneurons, motosynapses,
  • retinal neurons that feed visual information to the motor neurons : retina, sensorimotor_synapses

Top : Firing activity of the left (blue) and right(red) motor neurons

Bottom : retinal neuron firing activity while tracking a moving object during smooth pursuit

Strengths and Weaknesses

Strengths

  • Strong ability to implement biophysical detail in models
    • comparable level of detail to NEURON
  • easy to learn
  • able to implement a wide range of biological systems and phenomena
  • highly extensible

Drawbacks

  • poor support for parallelization means code can take much longer to run
  • knowledge of C++ in addition to python is required for most efficient code execution
  • for modelling convolutional spiking neural networks, packages with TensorFlow integration (Nengo) are better




Resources

Brian2 GitHub

Brian2 Documentation

ModelDB

  • code for Bono & Clopath(2017)

Comparison of Neural Network Simulators

NEURON Simulator

References


Bono, J., & Clopath, C. (2017). Modeling somatic and dendritic spike mediated plasticity at the single neuron and network level. Nature Communications, 8(1), 706–706. https://doi.org/10.1038/s41467-017-00740-z

Stimberg, M., Brette, R., & Goodman, D. F. M. (2019). Brian 2: an intuitive and efficient neural simulator. BioRxiv, 595710. https://doi.org/10.1101/595710

Stimberg, M., Goodman, D. F. M., Benichoux, V., & Brette, R. (2014). Equation-oriented specification of neural models for simulations. Frontiers in Neuroinformatics, 8, 6–6. https://doi.org/10.3389/fninf.2014.00006

Tikidji-Hamburyan, R. A., Narayana, V., Bozkus, Z., & El-Ghazawi, T. A. (2017). Software for Brain Network Simulations: A Comparative Study. Frontiers in Neuroinformatics, 11, 46–46. https://doi.org/10.3389/fninf.2017.00046