Creating Simulations in Python¶
The Pype9 package is organised into sub-packages loosely corresponding to each
pipeline (e.g. simulate
, plot
). The simulate
package contains the
sub-packages, neuron
and nest
, which provide the simulator-specific
calls to their respective backends.
All classes required to design and run simulations in these packages derive
from corresponding classes in the common
package, which defines a
consistent Public API across all backends. Therefore, code designed to
run on with one backend can be switched to another by simply changing the
package the simulator-specific classes are imported from (like in PyNN).
Note
The neuron
and nest
packages can be imported separately. Therefore,
only the simulator you plan to use needs to be available on your system.
Simulation Control¶
Simulation parameters such as time step, delay limits and seeds for pseudo random number generators are set within an instance of the Simulation class. Simulator objects (e.g. cells and connections) can only be instantiated within the context of an active Simulation instance, and there can only be one active Simulation instance at any time.
A Simulation is activated with the with
keyword
with Simulation(dt=0.1 * un.ms, seed=12345) as sim:
# Design simulation here
The simulation is advanced using the run
method of the Simulation
instance
with Simulation(dt=0.1 * un.ms, seed=12345) as sim:
# Create simulator objects here
sim.run(100.0 * un.ms)
this can be done in stages if states or parameters need to be updated mid-simulation
with Simulation(dt=0.1 * un.ms, seed=12345) as sim:
# Create simulator objects here
sim.run(50.0 * un.ms)
# Update simulator object parameters/state-variables
sim.run(50.0 * un.ms)
After the simulation context exits all objects in the simulator backend are destroyed (unless an exception is thrown) and only recordings can be reliably accessed from the “dead” Pype9 objects.
There are pseudo-random number generator (RNG) seeds that can be passed to
Simulation instances, seed
which is used to seed the RNGs that
generate random dynamic processes during the simulation (e.g. poisson
processes), and properties_seed
which is used to seed the RNG used to
generate probabilistic connectivity rules and the random distribution of cell
properties over populations. If seed
is None (the default) then it is
generated from the current timestamp, if properties_seed
is None then it is
derived from seed
.
Cell Simulations¶
NineML Dynamics classes can be translated into simulator cell objects using the
CellMetaClass class. A metaclass is class of classes, i.e. one whose
instantiation is itself a class, such as the type
class.
CellMetaClass instantiations derive from the Cell class and can
be used to represent different classes of neural models, such as Izhikevich or
Hodgkin-Huxley for example. From these Cell classes as many cell
instances (with their corresponding simulator objects) can be created as
required e.g:
# Create Izhikevich cell class by instantiating the CellMetaClass with a
# ninml.Dynamics Izhikevich model
Izhikevich = CellMetaClass('./izhikevich.xml#Izhikevich')
# Parameters and states of the cell class must be provided when the cells
# are instantiated.
# either as keyword args
izhi1 = Izhikevich(a=1, b=2, c=3, d=4, v=-65 * un.mV, u=14 * un.mV / un.ms)
# or from a nineml.DynamicsProperties object
izhi3 = Izhikevich('./izhikevich.xml#IzhikevichBurster')
If the specified Dynamics class has not been built before the CellMetaClass will automatically generate the required source code for the model, compile it, and load it into the simulator namespace. This can happen either inside or outside of an active Simulation instance. However, the cells objects themselves must be instantiated within a Simulation instance.
# The cell class can be created outside the simulation context
Izhikevich = CellMetaClass('./izhikevich.xml#Izhikevich')
with Simulation(dt=0.1 * un.ms) as sim:
# The cell object must be instantiated within the simulation context
izhi = Izhikevich(a=1, b=2, c=3, d=4, v=-65 * un.mV,
u=14 * un.mV / un.ms)
sim.run(1000.0 * un.ms)
The data can be recorded from every send port and state variable in the NineML
Dynamics class using the record
method of the Cell class. The
recorded data can then be accessed with the recording
method. The
recordings will be Neo format.
Izhikevich = CellMetaClass('./izhikevich.xml#Izhikevich')
with Simulation(dt=0.1 * un.ms) as sim:
izhi = Izhikevich(a=1, b=2, c=3, d=4, v=-65 * un.mV,
u=14 * un.mV / un.ms)
# Specify the variables to record
izhi.record('v')
sim.run(1000.0 * un.ms)
# Retrieve the recording
v = izhi.recording('v')
Transitions between regimes can be recorded using record_regime
and
retrieved using regime_epochs
Izhikevich = CellMetaClass('./izhikevich.xml#Izhikevich')
with Simulation(dt=0.1 * un.ms) as sim:
izhi = Izhikevich(a=1, b=2, c=3, d=4, v=-65 * un.mV,
u=14 * un.mV / un.ms)
# Specify the variables to record
izhi.record_regime()
sim.run(1000.0 * un.ms)
# Retrieve the regime changes
epochs = izhi.regime_epochs()
Data in Neo format can be “played” into receive ports of the Cell
neo_data = neo.PickleIO('./data/my_recording.neo.pkl').read()
Izhikevich = CellMetaClass('./izhikevich.xml#Izhikevich')
with Simulation(dt=0.1 * un.ms) as sim:
izhi = Izhikevich(a=1, b=2, c=3, d=4, v=-65 * un.mV,
u=14 * un.mV / un.ms)
# Play analog signal (must be of current dimension) into 'i_syn'
# analog-receive port.
izhi.play('i_syn', neo_data.analogsignals[0])
sim.run(1000.0 * un.ms)
States and parameters can be accessed and set using the attributes of the Cell objects
Izhikevich = CellMetaClass('./izhikevich.xml#Izhikevich')
with Simulation(dt=0.1 * un.ms) as sim:
izhi = Izhikevich(a=1, b=2, c=3, d=4)
sim.run(500.0 * un.ms)
# Update the membrane voltage after 500 ms to 20 mV
izhi.v = 20 * un.mV
sim.run(500.0 * un.ms)
Event ports can be connected between individual cells
Poisson = CellMetaClass('./poisson.xml#Poisson')
LIFAlphSyn = CellMetaClass('./liaf_alpha_syn.xml#LIFAlphaSyn')
with Simulation(dt=0.1 * un.ms) as sim:
poisson = Poisson(rate=10 * un.Hz, t_next=0.5 * un.ms)
lif = LIFAlphaSyn('./liaf_alpha_syn.xml#LIFAlphaSynProps')
# Connect 'spike_out' event-send port of the poisson cell to
# the 'spike_in' event-receive port on the leaky-integrate-and-fire
# cell
lif.connect(poisson, 'spike_out', 'spike_in')
sim.run(1000.0 * un.ms)
Network Simulations¶
Network simulations are specified in much the same way as Cell Simulations, with the exception that there is no metaclass for Networks (Network metaclasses will be added when the “Structure Layer” is introduced in NineML v2). Therefore, Network objects need to be instantiated within the simulation context.
with Simulation(dt=0.1 * un.ms) as sim:
# Network objects need to be instantiated within the simulation context
network = Network('./brunel/AI.xml#AI')
sim.run(1000.0 * un.ms)
During construction of the network, the NineML Populations and Projections are flattened into Component Array and Connection Group objects such that the synapse dynamics in the projection are included in the dynamics of the Component Array and each port connection is converted into a separate Connection Group of static connections.
To record data, the relevant component array needs to be accessed using the
component_array
or component_arrays
accessors of the network class.
Then as in the Cell Simulations case the record
method is used to
specify which variables to record and the recording
method is used to
access the recording after the simulation.
with Simulation(dt=0.1 * un.ms) as sim:
network = Network('./brunel/AI.xml#AI')
# 'spike_out' is explicitly connected in the connection so it is
# mapped to the global namespace of the flattened cell + synapses model
network.component_array('Exc').record('spike_out')
# State-variables of the cell dynamics are suffixed with '__cell'
network.component_array('Inh').record('v__cell')
# State-variables of synapses, in this case synapses from the
# 'Inhibition' projection, are prefixed with '__<projection-name>'
network.component_array('Exc').record('a__Inhibition')
sim.run(1000.0 * un.ms)
exc_spikes = network.component_array('Exc').recording('spike_out')
inh_v = network.component_array('Inh').recording('v__cell')
exc_inh_a = network.component_array('Exc').recording('a__Inhibition')
Note
During the cell and synapse flattening process the names of state variables
and unconnected ports will be suffixed with __cell
if they belong to the
population dynamics or __<my-projection>
if they belong to the synapse
of the a projection
Network models are simulated via integration with PyNN and therefore will run
on multiple processes using Open MPI (and Open MP_ for NEST) if the
calling Python script is run with mpirun
/mpiexec
.