Submodule: opexebo.analysis¶

Specialised neuroscience analysis functions

Overview¶

Positional analysis

`opexebo.analysis.spatial_occupancy`(time, …)	Generate an occpuancy map: how much time the animal spent in each location in the arena.
`opexebo.analysis.rate_map`(occupancy_map, …)	Calculate the spatially correlated firing rate map
`opexebo.analysis.rate_map_stats`(rate_map, …)	Calculate statistics of a rate map that depend on probability distribution function (PDF)
`opexebo.analysis.rate_map_coherence`(rate_map_raw)	Calculate coherence of a rate map
`opexebo.analysis.place_field`(firing_map, …)	Locate place fields on a firing map.
`opexebo.analysis.autocorrelation`(firing_map)	Calculate 2D spatial autocorrelation of a firing map.
`opexebo.analysis.grid_score`(aCorr, **kwargs)	Calculate gridness score for an autocorrelogram.
`opexebo.analysis.speed_score`(spike_times, …)	Calculate Speed score.

Angular analysis

`opexebo.analysis.angular_occupancy`(time, …)	Calculate angular occupancy from tracking angle and kwargs over (0,2*pi)
`opexebo.analysis.tuning_curve`(…)	Analogous to a RateMap - i.e.
`opexebo.analysis.tuning_curve_stats`(…)	Calculate statistics about a turning curve

Miscellanious

`opexebo.analysis.population_vector_correlation`(…)	Calculates the bin-wise correlation between two stacks of rate maps
`opexebo.analysis.theta_modulation_index`(…)	Calculate the 1-dimensional autocorrelation of the spike train and use it as a proxy score for how how reliably a cell fires at the same phase of the theta frequency.
`opexebo.analysis.calc_speed`(t, x[, y, …])	Calculate the speed of an animal given a list of times and positions

Experimental

`opexebo.analysis.border_coverage`(fields, …)	Calculate border coverage for detected fields.
`opexebo.analysis.border_score`(rate_map, …)	Calculate border score for a firing map

Detailed Documentation¶

opexebo.analysis.calc_speed(t, x, y=None, moving_average=None)¶

Calculate the speed of an animal given a list of times and positions

Speed is calculated element-wise, with an optional moving-average filter. Speed is not checked for non

Parameters:	t : np.ndarray 1d array containing time stamps x, y : np.ndarray 1d arrays containing positions. y is optional and may be excluded for 1d datasets moving_average: int, optional Apply a moving-average filter of this size to the speed data To exclude this filter, leave blank or set [None] (default)
Returns:	np.ndarray 1d array containing speed in units [position_unit / time_unit]

opexebo.analysis.spatial_occupancy(time, position, arena_size, **kwargs)¶

Generate an occpuancy map: how much time the animal spent in each location in the arena.

NOTES: This assumes that the positions have already been aligned and curated to remove NaNs. This is based on the expectation that it will primarily be used within the DataJoint framework, where the curation takes place at a much earlier stage.

Parameters:

time: np.ndarray

timestamps of position data

position: np.ndarray (x, [y])

1d or 2d array of positions at timestamps. If 2d, then row major, such that position[0] corresponds to all x; and position[1] to all y

arena_size: float or tuple of floats

Dimensions of arena (in cm)

For a linear track, length
For a circular arena, diameter
For a rectangular arena, length or (length, length)

arena_shape: {“square”, “rect”, “circle”, “line”}

Rectangular and square are equivalent. Elliptical or n!=4 polygons not currently supported. Defaults to Rectangular

bin_width: float

Bin size in cm. Default 2.5cm. If bin_width is supplied, limit must also be supplied. One of bin_width, bin_number, bin_edges must be provided

bin_number: int or tuple of int

Number of bins. If provided as a tuple, then (x_bins, y_bins). One of bin_width, bin_number, bin_edges must be provided

bin_edges: array-like

Edges of the bins. Provided either as edges or (x_edges, y_edges). One of bin_width, bin_number, bin_edges must be provided

limits: tuple or np.ndarray

(x_min, x_max) or (x_min, x_max, y_min, y_max) Provide concrete limits to the range over which the histogram searches Any observations outside these limits are discarded If no limits are provided, then use np.nanmin(data), np.nanmax(data) to generate default limits. As is standard in python, acceptable values include the lower bound and exclude the upper bound

debug: bool, optional

If true, print out debugging information throughout the function. Default False

Returns:

masked_map: np.ma.MaskedArray: Unsmoothed map of time the animal spent in each bin. Bins which the animal never visited are masked (i.e. the mask value is True at these locations)
coverage: float: Fraction of the bins that the animal visited. In range [0, 1]
bin_edges: ndarray or tuple of ndarray: x, or (x, y), where x, y are 1d np.ndarrays Here x, y correspond to the output histogram

opexebo.analysis.rate_map(occupancy_map, spikes_tracking, arena_size, **kwargs)¶

Calculate the spatially correlated firing rate map

The rate map is calculated by the number of spikes in a bin divided by the time the animal spent in that bin. Bins in which the animal spent no or very little time are masked such that the value is available for but typically excluded from future analyses.

The provided arena_size and bin_width **must* provide a number of bins such that the spike map has the same dimensions as the occupancy map.

Parameters:

occupancy_map: np.ndarray or np.ma.MaskedArray

Nx1 or Nx2 array of time spent by animal in each bin, with time in bins

spikes_tracking: np.ndarray

Nx2 or Nx3 array of spikes tracking: [t, x] or [t, x, y].

arena_size: float or tuple of floats

Dimensions of arena (in cm)

For a linear track, length
For a circular arena, diameter
For a rectangular arena, length or (length, length)

bin_width: float

Bin size in cm. Default 2.5cm. If bin_width is supplied, limit must also be supplied. One of bin_width, bin_number, bin_edges must be provided

bin_number: int or tuple of int

Number of bins. If provided as a tuple, then (x_bins, y_bins). One of bin_width, bin_number, bin_edges must be provided

bin_edges: array-like

Edges of the bins. Provided either as edges or (x_edges, y_edges). One of bin_width, bin_number, bin_edges must be provided

limits: tuple or np.ndarray

(x_min, x_max) or (x_min, x_max, y_min, y_max) Provide concrete limits to the range over which the histogram searches Any observations outside these limits are discarded If no limits are provided, then use np.nanmin(data), np.nanmax(data) to generate default limits. As is standard in python, acceptable values include the lower bound and exclude the upper bound

Returns:

rmap: np.ma.MaskedArray: 2D array, masked at locations of very low occupancy (t<1ms). Each cell gives the rate of neuron firing at that location.

See also

opexebo.general.accumulate_spatial
opexebo.general.bin_width_to_bin_number

Notes

BNT.+analyses.map()

opexebo.analysis.rate_map_stats(rate_map, time_map, debug=False)¶

Calculate statistics of a rate map that depend on probability distribution function (PDF)

Calculates information, sparsity and selectivity of a rate map. Calculations are done according to 1993 Skaggs et al. “An Information-Theoretic Approach to Deciphering the Hippocampal Code” paper. Another source of information is 1996 Skaggs et al. paper called “Theta phase precession in hippocampal neuronal populations and the compression of temporal sequences”.

Coherence is calculated based on RU Muller, JL Kubie “The firing of hippocampal place cells predicts the future position of freely moving rats”, Journal of Neuroscience, 1 December 1989, 9(12):4101-4110. The paper doesn’t provide information about how to deal with border values which do not have 8 well-defined neighbours. This function uses zero-padding technique.

Parameters:	rate_map: np.ma.MaskedArray Smoothed rate map: n x m array where cell value is the firing rate, masked at locations with low occupancy time_map: np.ma.MaskedArray time map: n x m array where the cell value is the time the animal spent in each cell, masked at locations with low occupancy Already smoothed
Returns:	rms: dict spatial_information_rate: float information rate [bits/sec] spatial_information_content: float spatial information content [bits/spike] sparsity: float see relevant literature (above) selectivity: float see relevant literature (above) peak_rate: float peak firing rate of smoothed map [Hz] mean_rate: float mean firing rate of smoothed map [Hz]

Notes

BNT.+analyses.mapStatsPDF(map)

BNT.+analyses.coherence(map)

opexebo.analysis.rate_map_coherence(rate_map_raw)¶

Calculate coherence of a rate map

Coherence is calculated based on RU Muller, JL Kubie “The firing of hippocampal place cells predicts the future position of freely moving rats”, Journal of Neuroscience, 1 December 1989, 9(12):4101-4110. The paper doesn’t provide information about how to deal with border values which do not have 8 well-defined neighbours. This function uses zero-padding technique.

Parameters:	rate_map_raw: np.ma.MaskedArray Non-smoothed rate map: n x m array where cell value is the firing rate, masked at locations with low occupancy
Returns:	coherence: float see relevant literature (above)

Notes

BNT.+analyses.coherence(map)

opexebo.analysis.grid_score(aCorr, **kwargs)¶

Calculate gridness score for an autocorrelogram.

Calculates a gridness score by expanding a circle around the centre field and calculating a correlation value of that circle with it’s rotated versions. The expansion is done up until the smallest side of the autocorrelogram. The function may also calculate grid statistics.

Gridness score value by itself is calculated as a maximum over a sliding mean of expanded circle. The width of the sliding window is given by a variable numGridnessRadii. This is done in order to keep gridness score the same as historical values (i.e. older versions of gridness score).

Other Parameters:
Parameters:	acorr: np.ndarray A 2D autocorrelogram.
Returns:	grid_score: float Always returns a gridness score value. It ranges from -2 to 2. 2 is more of a theoretical bound for a perfect grid. More practical value for a good grid is around 1.3. If function can not calculate a gridness score, NaN value is returned. grid_stats: dictionary grid_spacings: np.array Spacing of three adjacent fields closest to center in autocorr (in [bins] if keyword bin_width is not given) grid_spacing: float Nanmean of ‘spacings’ (in [bins] if keyword bin_width is not given) grid_orientations: np.array Orientation of three adjacent fields closest to center in autocorr (in [degrees]) grid_orientations_std: float Standard deviation of orientations % 60 grid_orientation: float Orientation of grid in [degrees] (mean of fields of 3 main axes) grid_positions: np.array [y,x] coordinates of six fields closest to center grid_ellipse: np.array Ellipse fit returning [x coordinate, y coordinate, major radius, minor radius, theta] grid_ellipse_aspect_ratio: float Ellipse aspect ratio (major radius / minor radius) grid_ellipse_thet: float Ellipse theta (corrected according to previous BNT standard) in [degrees]
	min_orientation: int See function “grid_score_stats” search_method: str Peak searching method, currently limited to either default or sep bin_width: float Size of the bins. Distance units will be returned in the same units If not provided, distance units will be retuned in units [bins]

See also

opexebo.analysis.placefield

Notes

BNT.+analyses.gridnessScore

opexebo.analysis.grid_score_stats(aCorr, mask, centre, **kwargs)¶

Calculate spatial characteristics of grid based on 2D autocorr

Parameters:

aCorr : np.array

2D Autocorrelation

mask : np.array

Mask (masked=True) of shape aCorr for masking center field and fringes of aCorr above best grid score radius

centre : np.array

Centre coordinate [y,x]

**kwargs :

min_orientation : int: Minimum difference in degrees that two neighbouring fields detected in 2D autocorrelation must have. If difference is below this threshold, discard the field that has larger distance from center
search_method : str: Peak searching method, currently limited to either default or sep
bin_width : float: Size of the bins. Distance units will be returned in the same units If not provided, distance units will be retuned in units [bins]

Returns:

grid_stats : dictionary

grid_spacings : np.array: Spacing of three adjacent fields closest to center in autocorr (in [bins])
grid_spacing : float: Nanmean of ‘spacings’ in [bins]
grid_orientations : np.array: Orientation of three adjacent fields closest to center in autocorr (in [degrees])
grid_orientations_std : float: Standard deviation of orientations % 60
grid_orientation : float: Orientation of grid in [degrees] (mean of fields of 3 main axes)
grid_positions : np.array: [y,x] coordinates of six fields closest to center
grid_ellipse : np.array: Ellipse fit returning [x coordinate, y coordinate, major radius, minor radius, theta]
grid_ellipse_aspect_ratio : float: Ellipse aspect ratio (major radius / minor radius)
grid_ellipse_theta : float: Ellipse theta (corrected according to previous BNT standard) in [degrees]

opexebo.analysis.autocorrelation(firing_map)¶

Calculate 2D spatial autocorrelation of a firing map.

Parameters:	firing_map: np.ndarray NxM matrix, smoothed firing map. map is not necessary a numpy array. May contain NaNs.
Returns:	acorr: np.ndarray Resulting correlation matrix, which is a 2D numpy array.

See also

opexebo.general.normxcorr2_general

Notes

BNT.+analyses.autocorrelation

opexebo.analysis.place_field(firing_map, **kwargs)¶

Locate place fields on a firing map.

Identifies place fields in 2D firing map. Placefields are identified by using an adaptive threshold. The idea is that we start with a peak value as the threshold. Then we gradually decrease the threshold until the field area doesn’t change any more or the area explodes (this means the threshold is too low).

Other Parameters:
Parameters:	firing_map: np.ndarray or np.ma.MaskedArray smoothed rate map. If supplied as an np.ndarray, it is assumed that the map takes values of np.nan at locations of zero occupancy. If supplied as an np.ma.MaskedArray, it is assumed that the map is masked at locations of zero occupancy
Returns:	fields: list of dict coords: np.ndarray Coordinates of all bins in the firing field peak_coords: np.ndarray Coordinates peak firing rate [y,x] centroid_coords: np.ndarray Coordinates of centroid (decimal) [y,x] area: int Number of bins in firing field. [bins] bbox: tuple Coordinates of bounding box including the firing field (y_min, x_min, y_max, y_max) mean_rate: float mean firing rate [Hz] peak_rate: float peak firing rate [Hz] map: np.ndarray Binary map of arena. Cells inside firing field have value 1, all other cells have value 0 fields_map : np.ndarray labelled integer image (i.e. background = 0, field1 = 1, field2 = 2, etc.)
	min_bins: int Fields containing fewer than this many bins will be discarded. Default 9 min_peak: float Fields with a peak firing rate lower than this absolute value will be discarded. Default 1 Hz min_mean: float Fields with a mean firing rate lower than this absolute value will be discarded. Default 0 Hz init_thresh: float Initial threshold to search for fields from. Must be in the range [0, 1]. Default 0.96 search_method: str Peak detection finding method. By default, use skimage.morphology.local_maxima Acceptable values are defined in opexebo.defaults. Not required if peak_coords are provided peak_coords: array-like List of peak co-ordinates to consider instead of auto detection. [y, x]. Default None
Raises:	ValueError Invalid input arguments NotImplementedError non-defined peack searching methods

Notes

BNT.+analyses.placefieldAdaptive

https://se.mathworks.com/help/images/understanding-morphological-reconstruction.html

opexebo.analysis.egocentric_occupancy(positions, head_angles, num_angles, num_bins, boundaries, arena_shape)¶

Generate a 2D histogram of (angle, border_distance) data . Can calculate either allocentric or egocentric values. To calculate allocentric (=fixed reference), set head_angles to None

Implementation is based on, but indepdent of, https://doi.org/10.1016/j.cub.2019.07.007

This is a moderately general function: it will work for either the tracking or spike data. However, it is sensible to enforce that the boundaries keyword uses the same value for both sides if the goal is to generate a meaningful ratemap

Parameters:

position: np.ndarray

array of (x, y) positions of animal

head_angle: np.ndarray -or- None

Angle [RADIANS] of the animal’s head with respect to the X axis Zero radians points along the x axis. Set to None to calculate the allocentric values

num_angles: int

Number of angles at which to calculate the wall distance

num_bins:

Number of bins to consider wall distance

boundaries: tuple

Some way of defining where the boundary of the arena is if arena_shape is square, then boundaries should be:

( (min_x, max_x), (min_y, max_y) )

elif arena_shape is circ, then boundaries should be:: ( (centre_x, centre_y), radius)

arena_shape: str

definition of the shape of the arena, ties in to format of boundaries Acceptable values are given in `opexebo.defaults.shapes_*

Returns:

hist: np.ndarray: 2D histogram of (wall_distance, angle) data hist[0, :] are all the distance bins for a single angle hist[:, 0] are all the angles for a single distance bin
distance_bins: np.ndarray: Distance bins, +1 element longer than hist[0, :].shape
angle_bins: np.ndarray: Angular bins, +1 element longer than hist[:, 0].shape

opexebo.analysis.angular_occupancy(time, angle, **kwargs)¶

Calculate angular occupancy from tracking angle and kwargs over (0,2*pi)

Parameters:	time : numpy.ndarray time stamps of angles in seconds angle : numpy array Head angle in radians Nx1 array bin_width : float, optional Width of histogram bin in degrees
Returns:	masked_histogram : numpy masked array Angular histogram, masked at angles at which the animal was never observed. A mask value of True means that the animal never occupied that angle. coverage : float Fraction of the bins that the animal visited. In range [0, 1] bin_edges : list-like x, or (x, y), where x, y are 1d np.ndarrays Here x, y correspond to the output histogram

Notes

This program is free software; you can redistribute it and/or modify it under the terms of the GNU General Public License as published by the Free Software Foundation; either version 3 of the License, or (at your option) any later version.

opexebo.analysis.tuning_curve(angular_occupancy, spike_angles, **kwargs)¶

Analogous to a RateMap - i.e. mapping spike activity to spatial position map spike rate as a function of angle

Parameters:

angular_occupancy : np.ma.MaskedArray

unsmoothed histogram of time spent at each angular range Nx1 array, covering the range [0, 2pi] radians Masked at angles of zero occupancy

spike_angles : np.ndarray

Mx1 array, where the m’th value is the angle of the animal (in radians) associated with the m’th spike

kwargs

bin_width : float: width of histogram bin in DEGREES Must match that used in calculating angular_occupancy In the case of a non-exact divisor of 360 deg, the bin size will be shrunk to yield an integer bin number.

Returns:

tuning_curve : np.ma.MaskedArray: unsmoothed array of firing rate as a function of angle Nx1 array

Notes

BNT.+analyses.turningcurve

This program is free software; you can redistribute it and/or modify it under the terms of the GNU General Public License as published by the Free Software Foundation; either version 3 of the License, or (at your option) any later version.

opexebo.analysis.tuning_curve_stats(tuning_curve, **kwargs)¶

Calculate statistics about a turning curve

STATUS : EXPERIMENTAL

Calculates various statistics for a turning curve. 1. Mean vector length of a head direction rate map. The value will range from 0 to 1. 0 means that there are so much dispersion that a mean angle cannot be described. 1 means that all data are concentrated at the same direction. Note that 0 does not necessarily indicate a uniform distribution. Calculation is based on Section 26.4, J.H Zar - Biostatistical Analysis 5th edition, see eq. 26.13, 26.14.

Parameters:

tuning_curve : np.ma.MaskedArray

Smoothed turning curve of firing rate as a function of angle Nx1 array

kwargs

percentile : float: Percentile value for the head direction arc calculation Arc is between two points with values around globalPeak * percentile. Value should be in range [0, 1]

Returns:

tcstat : dict

hd_score : float: Score for how strongly modulated by angle the cell is
hd_mvl : float: mean vector length
hd_peak_rate : float: Peak firing rate [Hz]
hd_mean_rate : float: Mean firing rate [Hz]
hd_peak_direction : float: Direction of peak firing rate [degrees]
hd_peak_direction_rad : float: Direction of peak firing rate
hd_mean_direction: float: Direction of mean firing rate [degrees]
hd_mean_direction_rad: float: Direction of mean firing rate
hd_stdev : float: Circular standard deviation [degrees]
halfCwInd : int: Indicies of at the start, end of the range defined by percentile (clockwise).
halfCcwInd : int: Indicies of at the start, end of the range defined by percentile (counter-clockwise).
halfCwRad : float: Angle of the start, end of the range defined by percentile
halfCcwRad : float: Angle of the start, end of the range defined by percentile
arc_angle_rad : float: Angle of the arc defined by percentile
arc_angle_rad : float: Angle of the arc defined by percentile

Notes

BNT.+analyses.tcStatistics

This program is free software; you can redistribute it and/or modify it under the terms of the GNU General Public License as published by the Free Software Foundation; either version 3 of the License, or (at your option) any later version.

opexebo.analysis.population_vector_correlation(stack_0, stack_1, **kwargs)¶

Calculates the bin-wise correlation between two stacks of rate maps

Each stack corresponds to a separate Task, or trial. Each layer is the ratemap for a single cell from that Task. The same units should be given in the same order in each stack.

Take a single column through the stack (i.e. 1 single bin/location in arena, with a firing rate for each cell), from each stack

In the original MatLab implementation, three output modes were supported

1D: (numYbins) - iterate over i
1. Take a 2D slice from each stack - all cells at all X positions at a
single Y position i
1. Reshape from 2D to 1D
2. Calculate the Pearson correlation coefficient between the two 1D
arrays
1. The value of pv_corr_1d[i] is the Pearson correlation coefficient
arising from Y position i
2D (numXbins x numYbins) - iterate over i
1. Take a 2D slice from each stack - all cells at all X positions at a
single Y position i
1. Calculate the 2D array (numXbins x numYbins) where the `[j,k]`th
value is the Pearson correlation coefficient between all observations at the j’th X location in stack_left and the k’th location in stack_right
1. The i’th row of pv_corr_2d is the DIAGONAL of the correlation matrix
i.e. where j==k i.e. the correlation of the the SAME location in each stack for all observations (numCells)
3D (numXbins x numYbins x iteration(=`numYbins`))

Same as 2D BUT take the whole correlation matrix, not the diagonal i.e. the full [j,k] correlatio between all X locations

A note on correlation in Numpy vs Matlab

Matlab’s corr(a, b) function returns the correlation of ab Numpy’s corrcoef function returns the normalised covariance matrix, which is:

aa ab ba aa

The normalised covariance matrix should be hermitian, but due to floating point accuracy, this is not actually guaranteed the MatLab function can be reproduced by taking either [0, 1] or [1,0] of the normalised covariance matrix.

If a, b are 2D matricies, then they should have shape (num_variables, num_observations) In the case of this function, where the iterator is over the Y values of the rate map, that means: (x_bins, num_cells)

Parameters:

stack_0: 3D array -or- list of 2D arrays

stack_1: 3D array -or- list of 2D arrays

stack_x[i] should return the i’th ratemap. This corresponds to a constructor like:

np.zeros(num_layers, y_bins, x_bins)

Alternatively, a list or tuple of 2D arrays may be supplied:: stack_x = (ratemap_0, ratemap_1, ratemap_2, …)

row_major: bool

Direction of iteration. If True, then each row is iterated over in turn and correlation is calculated per row. If False, then each column is iterated over in turn, and correlation is calculated per column. Default True (same behavior as in BNT)

Returns:

(p1, p2, p3)
p1: np.ndarray (1D, iterator x 1): Array of Pearson correlation coefficients. i’th value is given by the correlation of the i’th flattened slice of stack_0 to the i’th flattened slice of stack_1
p2: np.ndarray (2D, iterator x non-iterator): i’th row is the diagonal of the correlation matrix, i.e. the correlation of the same location (location i) in each stack, i.e. where j==k
p3: np.ndarray(3D, iterator x non-iterator x non-iterator): i’th array is the entire correlation matrix, rather than just the diagonal

Notes

BNT.+analyses.populationVectorCorrelation

opexebo.analysis.theta_modulation_index(spike_times, **kwargs)¶

Calculate the 1-dimensional autocorrelation of the spike train and use it as a proxy score for how how reliably a cell fires at the same phase of the theta frequency.

The measure is calculated as described in publication (with small modifications) Cacucci et al, Theta-Modulated Place-by-Direction Cellsin the Hippocampal Formation in the Rat, Journal of Neuroscience, 2004. doi: 10.1523/jneurosci.2635-04.2004

The calculation is very similar in practice, although somewhat less sophisicated, than the standard g(2) autocorrelation calculation.

For each spike, consider all spikes that follow within some of length P: Calculate the time delta to all spikes within P and bin as a histogram

Convert from number of occurrences to rate of occurrences (optional) Calculate Index as the contrast of the autocorrelation at Tau=(50,70ms) to Tau=(100,140ms) (assumed theta trough and peak respectively)

This indicator fails in two key areas:

It indicates nothing about the phase of theta, i.e. this is a scalar rather than vector calculation
It _assumes_ the value of theta for the animal in question, rather than calculating it directly

Parameters:	spike_times : np.ndarray Array of times at which spikes for a cell are detected, in [seconds]
Returns:	theta_modulation_index : float Contrast of autocorrelation hist : np.ndarray Histogram of time_delta to subsequent spikes within the next 0.5s bins : np.ndarray Bin edges of the histogram in [s]

opexebo.analysis.border_coverage(fields, arena_shape, **kwargs)¶

Calculate border coverage for detected fields.

STATUS : EXPERIMENTAL

This function calculates firing map border coverage that is further used in calculation of a border score. Coverage is calculated on a per-field basis, and is therefore dependent on the exact details of the place-field detection algorithm.

Additionally, please be aware of what top and bottom mean. This function treats top as the location where the y indicies are greatest. This is the standard convention in a GRAPH ((0,0) at bottom left, x increases to right and y increases to top). Images are conventionally plotted UPSIDE- DOWN, with y increasing towards the bottom.

Other Parameters:
Parameters:	fields: dict or list of dicts One dictionary per field. Each dictionary must contain the keyword field_map arena_shape: {“square”, “rect”, “circle”, “line”} Rectangular and square are equivalent. Elliptical or n!=4 polygons not currently supported.
Returns:	coverage: float Border coverage, ranges from 0 to 1.
	search_width: int rate_map and fields_map have masked values, which may occur within the region of border pixels. To mitigate this, we check rows/columns within search_width pixels of the border Default 8 walls: str or list of tuple Definition of which walls to consider for border coverage. Behaviour is different for different arena shapes * Rectangular arenas type str. The four walls are referred to as T, B, R, L for the standard cardinals. Case insensitive. You may include or exclude any Circular arenas Walls should be specified as a list of inclusive angular start-stop pairs. Angles given in degrees clockwise from top. Both 0 and 360 are permissible values. Start-stop may wrap through zero and may be overlapping with other pairs. Negative angles may be used to express counter-clockwise rotation, and will be automatically wrapped into the [0, 360] range. Example: [(0, 15), (180, 270), (315, 45), (30, 60)] Walls may ALTERNATIVELY be expressed as a string, equivalent to rectangular arenas. In this case, walls will be interpreted as follows T - (315, 45) B - (135, 225) L - (45, 135) R - (225, 315) In either case, defaults TRBL. Converting this (or other string) for circluar arrays is handled internally

See also

opexebo.analysis.placefield
opexebo.analysis.borderscore

Notes

BNT.+analyses.placefield
BNT.+analyses.borderScore
BNT.+analyses.borderCoverage

opexebo.analysis.border_score(rate_map, fields_map, fields, arena_shape, **kwargs)¶

Calculate border score for a firing map

STATUS : EXPERIMENTAL

TODO: Update to account for the use of Vadim’s placefield code TODO: Handle circular arenas (How?)

Calculates a border score for a firing rate map according to the article “Representation of Geometric Borders in the Entorhinal Cortex” by Solstad et. al. (Science 2008).

Border score is either -1 (no firing fields provided) or in the range [-1, 1]. The score reflects both the width of a field (what fraction of a single wall it touches), and the depth of a field (how far away from the field it extends) The highest score is returned for a field with maximum width and infinitesimal depth.

Consequently, a border score of 1 can only ever be returned given an infinite resolution. A perfect field in a typical 40x40 bin map has a maximum value of around 0.91. (limited by the width of the binning compared to the arena extent)

The score is only evaluated for the single firing field that has the greatest wall coverage. All other fields are ignored.

Parameters:

rate_map: np.ma.MaskedArray

rate map: N x M array where each cell value is the firing rate of the cell

fields_map: np.ma.MaskedArray

Integer array of labelled fields. Each cell value is a positive integer. Cells that are members of field have the value corresponding to field_id (non-zero positive unique integers. Not necessarily contiguous (e.g. 1, 2, 3, 5)). Cells that are not members of fields have value zero

fields: list of dict

List of dictionaries of firing fields. Each dictionary must, at least, contain the keyword “field_map”, yielding a binary map of that field within the overall arena

arena_shape: {“square”, “rect”, “circle”, “line”}

Rectangular and square are equivalent. Elliptical or n!=4 polygons not currently supported.

search_width: int, optional

rate_map and fields_map have masked values, which may occur within the region of border pixels. To mitigate this, we check rows/columns within search_width pixels of the border Default 8

walls: str or list of tuple, optional

Definition of which walls to consider for border coverage. Behaviour is different for different arena shapes * Rectangular arenas

type str. The four walls are referred to as T, B, R, L for the standard cardinals. Case insensitive. You may include or exclude any

Circular arenas Walls should be specified as a list of inclusive angular start-stop pairs. Angles given in degrees clockwise from top. Both 0 and 360 are permissible values. Start-stop may wrap through zero and may be overlapping with other pairs. Negative angles may be used to express counter-clockwise rotation, and will be automatically wrapped into the [0, 360] range. Example: [(0, 15), (180, 270), (315, 45), (30, 60)] Walls may ALTERNATIVELY be expressed as a string, equivalent to rectangular arenas. In this case, walls will be interpreted as follows
- T - (315, 45)
- B - (135, 225)
- L - (45, 135)
- R - (225, 315)

In either case, defaults TRBL. Converting this (or other string) for circluar arrays is handled internally

Returns:

score: float: Border score in the range [-1, 1]. A value of -1 is given when no fields are provided
coverage: float: Border coverage in the range [0, 1] The coverage is given for ONE border. Coverage for a specific wall is only available if a single wall is searched

See also

opexebo.analysis.placefield
opexebo.analysis.bordercoverage

Notes

BNT.+analyses.placefield
BNT.+analyses.borderScore
BNT.+analyses.borderCoverage

opexebo.analysis.speed_score(spike_times, tracking_times, tracking_speeds, **kwargs)¶

Calculate Speed score.

Speed score is a correlation between cell firing rate and animal speed. The Python version is based on BNT.+scripts.speedScore. At Edvard’s request, both the 2015 and 2016 scores are calculated. The primary difference is how the speed smoothing is implemented

Speed score originates in the following paper in Nature from Emiliano et al doi:10.1038/nature14622

The original Matlab script implemented - but as far as I can tell, did not (by default) use, a Kalman filter for smoothing the animal speed. Since it is not the default behaviour, I have not (yet) added that Kalman filter to opexebo. Its addition is contingent on the score similarity to BNT.

Summary:

The intention is to correlate (spike firing rate) with (animal speed)
Convert an N-length array of spike firing times into an M-length array of spike firing rates, where M is the same length as tracking times
Optional: smooth firing rates
Optional: smooth speeds
Optional: apply a bandpass filter to speeds
calculate the Pearson correlation coefficient between (speed), (firing rate)

Other Parameters:
Parameters:	spike_times: np.ndarray N-length array listing the times at which spikes occurred. [s] tracking_times: np.ndarray M-length array of time stamps of tracking frames tracking_speeds: np.ndarray M-length array of animal speeds at time stamps given in tracking_times
Returns:	scores: dict 2015: float 2016: float Variations on the speed score. ‘2015’ is based on the code in the paper above, but additionally including an upper speed filter ‘2016’ is a modification involving a slightly different approach to smoothing the firing rate data (lower_speed, upper_speed): list of floats Speed thresholds using in the bandpass filter Most useful in the case of the adaptive filter, because there is no other way to find out what was actually used.
	bandpass: str Type of bandpass filter applied to animal speeds. Acceptable values are “none” - No speed based filtering is applied “fixed” - a fixed lower and upper speed bound are used, based on keywords “lower_bound_speed”, “upper_bound_speed” “adaptive” - a fixed lower speed bound is used, based on keyword “lower_bound_speed”. An upper speed bound is determined based on keywords “upper_bound_time” and “speed_bandwidth” Default none lower_bound_speed: float Speed in [cm/s] used as the lower edge of the speed bandpass filter (“fixed” and “adaptive”). Default 2cm/s upper_bound_speed: float Speed in [cm/s] used as the upper edge of the speed bandpass filter (“fixed” only). Default 15cm/s upper_bound_time: float Duration in [s] used for determining the upper edge of the speed bandpass filter (“adaptive” only). Default 10s speed_bandwidth: float Range of speeds in [cm/s] used for determining the upper edge of the speed bandpass filter (“adaptive” only). Default 2cm/s sigma: float Standard deviation in [s] of Gaussian smoothing kernel for smoothing both speed and firing rate data. Default 0.5s debug: bool

Notes

BNT.+scripts.speedScore