Core functionality (dipde.internals)

External Population

class dipde.internals.externalpopulation.ExternalPopulation(firing_rate, record=False, **kwargs)[source]

Bases: object

External (i.e. background) source for connections to Internal Populations.

This class provides a background drive to internal population. It is used as the source argument to a connection, in order to provide background drive.

Parameters:

firing_rate : numeric, str

Output firing rate of the population. If numeric type, this defines a constant background generator; if str, it is interpreted as a SymPy function with independent variable ‘t’.

record : bool (default=False)

If True, a history of the output firing rate is recorded (firing_rate_record attribute).

**kwargs

Any additional keyword args are stored as metadata (metadata attribute).

Attributes

self.firing_rate_string (str) String representation of firing_rate input parameter.
self.metadata (dict) Dictonary of metadata, constructed by kwargs not parsed during construction.
self.firing_rate_record (list) List of firing rates recorded during Simulation.
self.t_record (list) List of times that firing rates were recorded during Simulation.
curr_firing_rate

Property that accesses the current firing rate of the population (Hz).

firing_rate(t)[source]

Firing rate of the population at time t (Hz).

initialize()[source]

Initialize the population at the beginning of a simulation.

Calling this method resets the recorder that tracks firing rate during a simulation. This method is called by the Simulation object ( initialization method), but can also be called by a user when defining an alternative time stepping loop.

initialize_firing_rate_recorder()[source]

Initialize recorder at the beginning of a simulation.

This method is typically called by the initialize method rather than on its own. It resets the lists that track the firing rate during execution of the simulation.

update()[source]

Update the population one time step.

This method is called by the Simulation object to update the population one time step. In turn, this method calls the update_firing_rate_recorder method to register the current firing rate with the recorder.

update_firing_rate_recorder()[source]

Record current time and firing rate, if record==True.

This method is called once per time step. If record is True, calling this method will append the current time and firing rate to the firing rate recorder.

Internal Population

class dipde.internals.internalpopulation.InternalPopulation(tau_m=0.02, v_min=-0.1, v_max=0.02, dv=0.0001, record=True, curr_firing_rate=0.0, update_method='approx', approx_order=None, tol=1e-12, norm=inf, **kwargs)[source]

Bases: object

Population density class

This class encapulates all the details necessary to propagate a population density equation driven by a combination of recurrent and background connections. The voltage (spatial) domain discretization is defined by linear binning from v_min to v_max, in steps of dv (All in units of volts). The probability densities on this grid are recorded pv, and must always sum to 1.

Parameters:

tau_m : float (default=.02)

Time constant (unit: 1/sec) of neuronal population.

v_min : float (default=-.1)

Minimum of voltage domain (unit: volt).

v_max : float (default=.02)

Maximum of voltage domain (Absorbing boundary), i.e spiking threshold (unit: volt).

dv : float (default=.0001)

Voltage domain discritization size (unit: volt).

record : bool (default=False)

If True, a history of the output firing rate is recorded (firing_rate_record attribute).

curr_firing_rate : float (default=0.0

Initial/Current firing rate of the population (unit: Hz).

update_method : str ‘approx’ or ‘exact’ (default=’approx’)

Method to update pv (exact can be quite slow).

approx_order : int or None (default=None)

Maximum Taylor series expansion order to use when computing update to pv.

tol : float (default=1e-12)

Error tolerance used when computing update to pv.

norm : non-zero int, np.inf, -np.inf, or ‘fro’ (default=np.inf)

Vector norm used in computation of tol.

**kwargs

Any additional keyword args are stored as metadata (metadata attribute).

Attributes

self.edges (np.array) Vector defining the boundaries of voltage bins.
self.pv (np.array) Vector defining the probability mass in each voltage bin (self.pv.sum() = 1).
self.firing_rate_record (list) List of firing rates recorded during Simulation.
self.t_record (list) List of times that firing rates were recorded during Simulation.
self.leak_flux_matrix (np.array) Matrix that defines the flux between voltage bins.
get_total_flux_matrix()[source]

Create a total flux matrix by summing presynaptic inputs and the leak matrix.

initialize()[source]

Initialize the population at the beginning of a simulation.

In turn, this method:

  1. Initializes the voltage edges (self.edges) and probability mass in each bin (self.pv),
  2. Creates an initial dictionary of inputs into the population, and
  3. Resets the recorder that tracks firing rate during a simulation.

This method is called by the Simulation object (initialization method), but can also be called by a user when defining an alternative time stepping loop.

initialize_edges()[source]

Initialize self.edges and self.leak_flux_matrix attributes.

This method initializes the self.edges attribute based on the v_min, v_max, and dv settings, and creates a corresponding leak flux matrix based on this voltage discretization.

initialize_firing_rate_recorder()[source]

Initialize recorder at the beginning of a simulation.

This method is typically called by the initialize method rather than on its own. It resets the lists that track the firing rate during execution of the simulation.

initialize_probability()[source]

Initialize self.pv to delta-distribution at v=0.

initialize_total_input_dict()[source]

Initialize dictionary of presynaptic inputs at beginning of simulation

This method is typically called by the initialize method rather than on its own. It creates a dictionary of synaptic inputs to the population.

n_bins

Number of probability mass bins.

n_edges

Number of probability mass bin edges.

plot_probability_distribution(ax=None)[source]

Convenience method to plot voltage distribution.

Parameters:

ax : None or matplotlib.pyplot.axes object (default=None)

Axes object to plot distribution on. If None specified, a figure and axes object will be created.

source_connection_list

List of all connections that are a source for this population.

update()[source]

Update the population one time step.

This method is called by the Simulation object to update the population one time step. In turn, this method:

  1. Calls the update_total_input_dict method to gather the current strengths of presynaptic input populations,
  2. Calls the update_propability_mass method to propagate self.pv one time-step,
  3. Calls the update_firing_rate method to compute the firing rate of the population based on flux over threshold, and
  4. Calls the update_firing_rate_recorder method to register the current firing rate with the recorder.
update_firing_rate()[source]

Update curr_firing_rate attribute based on the total flux of probability mass across threshold.

update_firing_rate_recorder()[source]

Record current time and firing rate, if record==True.

This method is called once per time step. If record is True, calling this method will append the current time and firing rate to the firing rate recorder.

update_propability_mass()[source]

Create a total flux matrix, and propogate self.pv one time-step.

update_total_input_dict()[source]

Update the input dictionary based on the current firing rates of presynaptic populations.

Connection

Module containing Connection class, connections between source and target populations.

class dipde.internals.connection.Connection(source, target, nsyn, **kwargs)[source]

Bases: object

Class that connects dipde source population to dipde target population.

The Connection class handles all of the details of propogating connection information between a source population (dipde ExtermalPopulation or InternalPopulation) and a target population (dipde InternalPopulation). Handles delay information via a delay queue that is rolled on each timestep, and reuses connection information by using a ConnectionDistribution object with a specific signature to avoid duplication among identical connections.

Parameters:

source : InternalPopulation or ExternalPopulation

Source population for connection.

target : InternalPopulation

Target population for connection.

nsyn : int

In-degree of connectivity from source to target.

weights : list

Weights defining synaptic distribution (np.ndarray).

probs : list (same length as weights, and sums to 1)

Probabilities corresponding to weights.

delay: float (default=0)

Transmission delay (units: sec).

metadata: Connection metadata, all other kwargs

curr_delayed_firing_rate

Current firing rate of the source (float).

Property that accesses the firing rate at the top of the delay queue, from the source population.

initialize()[source]

Initialize the connection at the beginning of a simulation.

Calling this method:

  1. Initializes a delay queue used to store values of inputs in a last-in-first-out rolling queue.
  2. Creates a connection_distribution object for the connection, if a suitable object is not already registered with the simulation-level connection distribution collection.

This method is called by the Simulation object (initialization method), but can also be called by a user when defining an alternative time stepping loop.

initialize_connection_distribution()[source]

Create connection distribution, if necessary.

If the signature of this connection is already registered to the simulation-level connection distribution collection, it is associated with self. If not, it adds the connection distribution to the collection, and associates it with self.

initialize_delay_queue()[source]

Initialiaze a delay queue for the connection.

The delay attribute of the connection defines the transmission delay of the signal from the souce to the target. Firing rate values from the source population are held in a queue, discretized by the timestep, that is rolled over once per timestep. if source is an ExternalPopulation, the queue is initialized to the firing rate at t=0; if the source is an InternalPopulation, the queue is initialized to zero.

update()[source]

Update Connection, called once per timestep.

Connection Distribution

Module containing ConnectionDistribution class, organizing connection details.

class dipde.internals.connectiondistribution.ConnectionDistribution(edges, weights, probs, sparse=True)[source]

Bases: object

Container for matrix that computes element flux, and associated methods.

A ConnectionDistribution object contains a matrix used to compute the time propogation of the voltage probability distribution, and a vector used to compute the flux over threshold at the end of the timestep. It is unique up to the descretization of the target voltage distribution, and the specific synaptic weight distribution of the connection. Therefore, whenever possible, a ConnectionDistribution object is reused for multiple connections with identical signatures, to reduce memory consumption.

Parameters:

edges: np.ndarray

Voltage bin discretization of target population.

weights: np.ndarray

Weights defining the discrete synaptic distribution.

probs : np.ndarray

Probabilities corresponding to weights.

Attributes

self.threshold_flux_vector (np.ndarray) Vector used to compute over-threshold flux.
self.flux_matrix (np.ndarray) Matrix used to propagate voltage distribution.
initialize()[source]

Initialize connection distribution.

Initialization creates the flux_matrix and threshold_flux_vector. Implemented lazily via a try catch when flux matrix is requested, that does not exist.

signature

A unique key used to organize ConnectionDistributions.

Connection Distribution Collection

class dipde.internals.connectiondistributioncollection.ConnectionDistributionCollection[source]

Bases: dict

Container that organizes connection components for a simulation, to reduce redundancy.

In a simulation, connections that share the same weights and probabilities, as well as the same target bin edges, can make use of the same flux_matrix and threshold_flux_vector. This can significantly improve the overall memory efficiency of the simulation. To facilitate this, each simulation creates a ConnectionDistributionCollection object that indexes the ConnectionDistribution objects according to their signature, and re-uses the for multiple connections.

add_connection_distribution(cd)[source]

Try and add a ConnectionDistribution object, if signature not already used.

Simulation

class dipde.internals.simulation.Simulation(population_list, connection_list, verbose=True)[source]

Bases: object

Initialize and run a dipde simulation.

The Simulation class handles the initialization of population and connection objects, and provides a convenience time stepping loop to drive a network simulation. Typical usage involves the create of populations and connections, construction of a simulation object, and then call to simulation.run()

Parameters:

population_list : list of ExternalPopulation, InternalPopulation objects

List of populations to include in simulation.

connection_list : list of Connection objects

List of connections to include in simulation.

verbose : bool (default=True)

Setting True prints current time-step at each update evaluation.

initialize(t0=0.0)[source]

Initialize simulation, populations, and connections.

This function is typically called by the self.run() method, however can be called independently if defining a new time stepping loop.

Parameters:

t0 : float (default=0.)

Simulation start time (unit=seconds).

run(t0=0.0, dt=0.001, tf=0.1)[source]

Main iteration control loop for simulation

The time step selection must be approximately of the same order as dv for the internal populations, if the ‘approx’ time stepping method is selected.

Parameters:

t0 : float (default=0.)

Simulation start time (unit=seconds), passed to initialize call.

tf : float (default=.1)

Simulation end time (unit=seconds).

dt : float (default=0.001)

Time step (unit=seconds).

Utilities

dipde.internals.utilities.approx_update_method_order(J, pv, dt=0.0001, approx_order=2)[source]

Approximate the effect of a matrix exponential, truncating Taylor series at order ‘approx_order’.

dipde.internals.utilities.approx_update_method_tol(J, pv, tol=2.2e-16, dt=0.0001, norm='inf')[source]

Approximate the effect of a matrix exponential, with residual smaller than tol.

dipde.internals.utilities.assert_probability_mass_conserved(pv)[source]

Assert that probability mass in control nodes sums to 1.

dipde.internals.utilities.descretize(distribution, N, loc=0, scale=1, shape=[])[source]

Compute a discrete approzimation to a scipy.stats continuous distribution.

dipde.internals.utilities.exact_update_method(J, pv, dt=0.0001)[source]

Given a flux matrix, pdate probabilty vector by computing the matrix exponential.

dipde.internals.utilities.flux_matrix(v, w, lam, p=1)[source]

Compute a flux matrix for voltage bins v, weight w, firing rate lam, and probability p.

dipde.internals.utilities.fraction_overlap(a1, a2, b1, b2)[source]

Calculate the fractional overlap between range (a1,a2) and (b1,b2).

Used to compute a reallocation of probability mass from one set of bins to another, assuming linear interpolation.

dipde.internals.utilities.get_v_edges(v_min, v_max, dv)[source]

Construct voltage edegs, linear discretization.

dipde.internals.utilities.get_zero_bin_list(v)[source]

Return a list of indices that map flux back to the zero bin.

dipde.internals.utilities.leak_matrix(v, tau)[source]

Given a list of edges, construct a leak matrix with time constant tau.

dipde.internals.utilities.redistribute_probability_mass(A, B)[source]

Takes two ‘edge’ vectors and returns a 2D matrix mapping each ‘bin’ in B to overlapping bins in A. Assumes that A and B contain monotonically increasing edge values.