**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.

Brian2 is simple to install

```
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`

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
```

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 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

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

- here, synaptic weight

`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

- 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

Neurons spike when their membrane voltage

**($v$)**changes enough to reach a certain thresholdThis 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$)**

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**

- 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} $$
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 variablein the complex version,

`volt`

indicates that $v$ is now defined with units,`I : amp`

defines an state variable`I`

with the unit amps, and`tau : 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

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

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

- 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]:

**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]:

- code for this model can be found on ModelDB

- 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

- 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

- 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

- 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

- code for Bono & Clopath(2017)

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