Report

BRIAN Simulator 11/4/11 NEURON is cool, but… • …it’s not suited particularly well for large network simulations • What if you want to look at properties of thousands of neurons interacting with one another? • What about changing properties of synapses through time? BRIAN Simulator • BRIAN allows for efficient simulations of large neural networks • Includes nice routines for setting up random connectivity, recording spike times, changing synaptic weights as a function of activity • www.briansimulator.org – Should be able to run from the unzipped brian directory – http://www.briansimulator.org/docs/installation.html Make sure it works! from brian import * brian_sample_run() Building blocks of BRIAN • Model consists of: – Equations – NeuronGroup – a group of neurons which obeys those equations – Connection – a way to define connection matrices between neurons – SpikeMonitor (and others) – a way to measure properties of the simulated network My First BRIAN Model • To start, let’s walk through the example code on briansimulator.org homepage • Copy this text into IDLE (or whichever .py editor you’re using) • Make sure it runs first! from brian import * eqs = ''' dv/dt = (ge+gi-(v+49*mV))/(20*ms) : volt dge/dt = -ge/(5*ms) : volt dgi/dt = -gi/(10*ms) : volt ''' P = NeuronGroup(4000, eqs, threshold=50*mV, reset=-60*mV) P.v = -60*mV Pe = P.subgroup(3200) Pi = P.subgroup(800) Ce = Connection(Pe, P, 'ge', weight=1.62*mV, sparseness=0.02) Ci = Connection(Pi, P, 'gi', weight=9*mV, sparseness=0.02) M = SpikeMonitor(P) run(1*second) raster_plot(M) show() My First BRIAN Model • • • First – we set up the model Notice the units! These are important in BRIAN, they help ensure everything you’re modeling makes sense Equations are written as strings, these are executed by the differential equation solver in BRIAN – Unit names/abbreviations are reserved keywords in BRIAN • Create our NeuronGroup using this model – Define number of neurons, model used, threshold & reversal potentials • Questions: – What is the reversal potential here? – How does this model differ from HH? from brian import * eqs = ''' dv/dt = (ge+gi-(v+49*mV))/(20*ms) : volt dge/dt = -ge/(5*ms) : volt dgi/dt = -gi/(10*ms) : volt ''' P = NeuronGroup(4000, eqs, threshold=50*mV, reset=-60*mV) P.v = -60*mV Pe = P.subgroup(3200) Pi = P.subgroup(800) Ce = Connection(Pe, P, 'ge', weight=1.62*mV, sparseness=0.02) Ci = Connection(Pi, P, 'gi', weight=9*mV, sparseness=0.02) M = SpikeMonitor(P) run(1*second) raster_plot(M) show() My First BRIAN Model • Initialize the neurons to starting potentials • Connect them – Here, constant weight, random connectivity from each subpopulation (excitatory & inhibitory) to all neurons • Which state variable to propagate to in “target” neuron: ‘ge’ for excitatory synapses, ‘gi’ for inhibitory from brian import * eqs = ''' dv/dt = (ge+gi-(v+49*mV))/(20*ms) : volt dge/dt = -ge/(5*ms) : volt dgi/dt = -gi/(10*ms) : volt ''' P = NeuronGroup(4000, eqs, threshold=50*mV, reset=-60*mV) P.v = -60*mV Pe = P.subgroup(3200) Pi = P.subgroup(800) Ce = Connection(Pe, P, 'ge', weight=1.62*mV, sparseness=0.02) Ci = Connection(Pi, P, 'gi', weight=9*mV, sparseness=0.02) M = SpikeMonitor(P) run(1*second) raster_plot(M) show() My First BRIAN Model • Hook up your electrophysiology equipment (here, measure spike times) – Can also record population rate, ISI – search documentation for Monitor – Only this info is saved from the simulation • Run the simulation! • Plot – Other plotting tools (hist_plot, Autocorrelograms, pylab, etc) from brian import * eqs = ''' dv/dt = (ge+gi-(v+49*mV))/(20*ms) : volt dge/dt = -ge/(5*ms) : volt dgi/dt = -gi/(10*ms) : volt ''' P = NeuronGroup(4000, eqs, threshold=50*mV, reset=-60*mV) P.v = -60*mV Pe = P.subgroup(3200) Pi = P.subgroup(800) Ce = Connection(Pe, P, 'ge', weight=1.62*mV, sparseness=0.02) Ci = Connection(Pi, P, 'gi', weight=9*mV, sparseness=0.02) M = SpikeMonitor(P) run(1*second) raster_plot(M) show() My First BRIAN Model from brian import * eqs = ''' dv/dt = (ge+gi-(v+49*mV))/(20*ms) : volt dge/dt = -ge/(5*ms) : volt dgi/dt = -gi/(10*ms) : volt ''' P = NeuronGroup(4000, eqs, threshold=-50*mV, reset=60*mV) P.v = -60*mV Pe = P.subgroup(3200) Pi = P.subgroup(800) Ce = Connection(Pe, P, 'ge', weight=1.62*mV, sparseness=0.02) Ci = Connection(Pi, P, 'gi', weight=-9*mV, sparseness=0.02) M = SpikeMonitor(P) R = PopulationRateMonitor(P,.01*second) H = ISIHistogramMonitor(P,bins = [0*ms,5*ms,10*ms,15*ms,20*ms,25*ms,30*ms,35*ms,40*m s]) run(1*second) raster_plot(M) hist_plot(H) plt.figure() plt.plot(R.times,R.rate) plt.xlabel(‘time (s)’) plt.ylabel(‘firing rate (Hz)’) show() Exercise: Can you make HW 4’s networks in BRIAN? • Try with both HH-style and integrate and fire models • Use documentation at briansimulator.org for help (especially Connection and Equations) – Hint: search HodgkinHuxley and Library Models • Let’s use simplified exponential synapses – Time constant tau = 3 ms for excitatory, 6 ms for inhibitory Feedforward Inhibition Feedback Inhibition Exercise: Can you make HW 4’s networks in BRIAN? Feedforward Inhibition Feedback Inhibition Spike timing dependent plasticity Final weight time Synapse number Count • Sets up network of Poisson spiking neurons, all projecting to a single neuron • Each synapse implements an STDP rule • By end of simulation, synapses become either strong or weak Firing rate (Hz) stdp1.py Final weight Final weight time Synapse number Count • Only difference is that we’ve changed how pre-before-post and post-before-pre synapses are weighted (and, in model code, stdp2 uses ExponentialSTDP, to make things easier) Firing rate (Hz) stdp2.py Final weight