SpikeSliceStack

The SpikeSliceStack class stores event-aligned slices of raw spike train data. It holds a list of SpikeData objects (one per trial or event window) and provides methods for converting to raster arrays of shape (num_units, num_bins, num_slices) and computing trial-averaged statistics.

class spikelab.spikedata.spikeslicestack.SpikeSliceStack(data_obj=None, times_start_to_end=None, time_peaks=None, time_bounds=None, spike_stack=None, neuron_attributes=None, drop_slice_attributes=True)

Bases: object

A list of SpikeData objects, one per slice, with spike-based comparison capabilities.

U is units (neurons) and S is slices (bursts, events, etc). Construct from either a single SpikeData with time specifications, or directly from a pre-built list of SpikeData objects.

Parameters:

• data_obj (SpikeData or None) – A SpikeData object to slice. Provide either this or spike_stack, not both.

• times_start_to_end (list or None) – Each entry is a tuple (start, end) representing the start and end times of a desired slice. Each tuple must have the same duration.

• time_peaks (list or None) – List of times as int or float where there is a burst peak or stimulation event. Must be paired with time_bounds. Alternative to times_start_to_end.

• time_bounds (tuple or None) – Single tuple (left_bound, right_bound). For example, (250, 500) means 250 ms before peak and 500 ms after peak. Must be paired with time_peaks.

• spike_stack (list or None) – List of SpikeData objects, one per slice. All must have the same number of units. Spike times must be relative to the slice (0-based or event-centered via start_time), not absolute recording times. Provide either this or data_obj, not both.

• neuron_attributes (list or None) – List of attribute dicts, one per unit. If None, inherited from data_obj when available.

• drop_slice_attributes (bool) – If True (default), neuron_attributes are removed from individual SpikeData slices after construction. The shared copy is stored at neuron_attributes. This avoids duplicating large per-unit data (e.g. waveform templates) across every slice. Set to False to keep per-slice attributes.

spike_stack

List of SpikeData objects, one per slice. Spike times are relative to the slice window. For 0-based slices, times run from 0 to duration. For event-centered slices, times run from -pre_ms to +post_ms with t=0 at the event. Use self.times for absolute recording time positions.

Type:

list

times

List of (start, end) time bounds for each slice in absolute recording time, sorted chronologically. Length equals S. Example: [(100, 350), (500, 750), (1000, 1250)].

Type:

list

N

Number of units.

Type:

int

neuron_attributes

List of attribute dicts, one per unit. None if not provided. When drop_slice_attributes is True (default), this is the only copy and individual slices will have neuron_attributes set to None.

Type:

list or None

__init__(data_obj=None, times_start_to_end=None, time_peaks=None, time_bounds=None, spike_stack=None, neuron_attributes=None, drop_slice_attributes=True)

subslice(slices)

Extract a subset of slices from the spike stack.

Parameters:

slices (int or list) – Slice index or list of slice indices to extract.

Returns:

New SpikeSliceStack containing only the

specified slices. Shape changes from S to S_trimmed. All units and neuron_attributes are carried over.

Return type:

result (SpikeSliceStack)

subset(units, by=None)

Extract a subset of units from every slice in the spike stack.

Parameters:

• units (int, str, or list) – Unit indices to extract. If by is None, must be int(s). If by is set, values to match in neuron_attributes.

• by (str or None) – If set, select units by this neuron_attribute key instead of by index.

Returns:

New SpikeSliceStack containing only the

specified units across all slices. All slices and neuron_attributes are carried over.

Return type:

result (SpikeSliceStack)

Notes

• Units are included in the output in the order they appear in the train (ascending index order), not the order listed in units.

• If IDs are not unique (when using by), every matching neuron is included.

subtime_by_index(start_idx, end_idx)

Trim each slice to a sub-window specified by millisecond indices.

Indices are measured from the start of each slice (1 index = 1 ms). Trims along the time axis while preserving all slices and units.

Parameters:

• start_idx (int) – Start index in ms from slice start (inclusive). Supports negative indexing.

• end_idx (int) – End index in ms from slice start (exclusive). Supports negative indexing.

Returns:

New SpikeSliceStack where each slice is

trimmed to the corresponding absolute time window. Absolute spike times are preserved (not shifted). self.times is updated to reflect the new absolute time bounds.

Return type:

result (SpikeSliceStack)

Notes

• Indices are relative to each slice’s own start (index 0 = slice start ms). They are converted to absolute recording times internally before trimming.

• Original absolute timestamps are preserved. To get shifted-to-zero timestamps, create a new SpikeSliceStack.

• All slices and neuron_attributes are carried over from the original.

to_raster_array(bin_size=1.0, absolute_times=False)

Convert the spike stack into a 3D raster array of shape (N, T, S).

Each slice is rasterized with the given bin size, producing a spike count matrix where entry (n, t, s) is the number of spikes unit n fired in time bin t of slice s.

Parameters:

• bin_size (float) – Time bin size in ms (default 1.0).

• absolute_times (bool) – If False (default), time bin 0 corresponds to the start of each slice (0-based). If True, each slice’s spikes are offset by its absolute start time from self.times, so bin indices reflect the original recording position. The T dimension is sized to cover the full time span from the earliest slice start to the latest slice end. Caution: this can produce very large arrays when the recording span is long and bin_size is small.

Returns:

3D array of shape (N, T, S) with

non-negative integer spike counts. When absolute_times is True, T covers the full recording span and all slices share the same time axis.

Return type:

raster_stack (np.ndarray)

compute_frac_active(min_spikes=2)

Compute the fraction of slices each unit is active in.

A unit counts as active in a slice if it has at least min_spikes spikes within that slice’s time window.

Parameters:

min_spikes (int) – Minimum number of spikes for a unit to count as active in a slice (default: 2).

Returns:

1-D array of shape (U,) with the

fraction of slices each unit is active in (values in [0, 1]).

Return type:

frac_active (np.ndarray)

Notes

• The returned array can be passed as frac_active to RateSliceStack.order_units_across_slices, RateSliceStack.get_slice_to_slice_unit_corr_from_stack, SpikeSliceStack.order_units_across_slices, or SpikeSliceStack.get_slice_to_slice_unit_comparison to override their internal activity calculation.

• SpikeData.get_frac_active produces a compatible (U,) array based on burst edges and can be used in the same way.

order_units_across_slices(agg_func='median', timing='median', min_spikes=2, min_frac_active=0.0, frac_active=None, timing_matrix=None)

Reorder units by their typical spike timing across slices.

For each unit in each slice, computes a representative spike time (median, mean, or first spike) relative to the slice’s time origin. These per-slice values are aggregated across slices to obtain a single typical timing per unit. Units are then sorted by this value from earliest to latest and optionally split into a highly-active group and a low-activity group.

Parameters:

• agg_func (str) – How to aggregate per-slice timing values across slices. "median" (default) or "mean".

• timing (str) – Which spike time to extract per unit per slice. "median" — median spike time within the slice (default). "mean" — mean spike time within the slice. "first" — first spike time (onset latency). Ignored when timing_matrix is provided.

• min_spikes (int) – Minimum number of spikes for a unit to count as active in a slice (default: 2). Ignored when timing_matrix is provided.

• min_frac_active (float or None) – Minimum fraction of slices a unit must be active in to be placed in the highly-active group. 0.0 or None (default: 0.0) skips the split entirely and places all units in the highly-active group without computing activity fractions.

• frac_active (np.ndarray or None) – Optional pre-computed fraction-active array of shape (U,) to override the internal calculation for the group split. Only used when min_frac_active > 0. Compatible sources: SpikeSliceStack.compute_frac_active and SpikeData.get_frac_active (frac_per_unit output).

• timing_matrix (np.ndarray or None) – Optional pre-computed (U, S) timing matrix from get_unit_timing_per_slice. When provided, timing and min_spikes are ignored and this matrix is used directly.

Returns:

Two SpikeSliceStack objects

(highly_active, low_active) with units reordered by typical timing. The low-activity stack is None when the group is empty.

unit_ids_in_order (tuple): Two ndarray

(highly_active, low_active) of original unit indices in the reordered sequence.

unit_std (tuple): Two ndarray (highly_active, low_active)

of standard deviation of per-slice timing values. Lower values indicate more consistent timing across slices.

unit_peak_times_ms (tuple): Two ndarray

(highly_active, low_active) of the aggregated typical timing in milliseconds relative to slice start. NaN for units with no active slices.

unit_frac_active (tuple): Two ndarray

(highly_active, low_active) of the fraction of slices each unit was active in.

Return type:

reordered_stacks (tuple)

Notes

• Call get_unit_timing_per_slice first to pre-compute the timing matrix if you want to reuse it across multiple calls (e.g. rank_order_correlation and this method).

• When frac_active is None and min_frac_active > 0, activity fraction is computed via compute_frac_active.

• Analogous to RateSliceStack.order_units_across_slices but operates on raw spike trains instead of firing rate curves.

apply(func, *args, **kwargs)

Apply a function to each SpikeData in the stack and return stacked results.

Calls func(sd, *args, **kwargs) on every slice and stacks the outputs into a single numpy array with a new leading axis of size S (number of slices).

Parameters:

• func (callable) – Function that accepts a SpikeData as its first argument and returns a numeric value (scalar, 1-D, or 2-D array). Output shape must be consistent across all slices.

• *args – Additional positional arguments forwarded to func.

• **kwargs – Additional keyword arguments forwarded to func.

Returns:

Stacked results with shape (S, ...).

Return type:

result (np.ndarray)

Notes

• Intended for use with stacks built by SpikeData.frames, SpikeData.align_to_events, SpikeData.spike_shuffle_stack, or SpikeData.subset_stack. Pair with shuffle_z_score, shuffle_percentile, slice_trend, or slice_stability from utils to interpret the results.

unit_to_unit_comparison(metric='ccg', delt=20.0, bin_size=1.0, max_lag=350, n_jobs=-1)

Compute pairwise unit-to-unit similarity within each slice using spike-based metrics.

For each slice, computes a (U, U) similarity matrix between all unit pairs, then stacks the results into a PairwiseCompMatrixStack (U, U, S).

Parameters:

• metric (str) – Similarity metric to use. "ccg" for cross-correlogram on binned rasters (default), "sttc" for spike time tiling coefficient.

• delt (float) – STTC time window in milliseconds (default: 20.0). Only used when metric is "sttc".

• bin_size (float) – Bin size in milliseconds for the binary raster (default: 1.0). Only used when metric is "ccg".

• max_lag (float) – Maximum lag in milliseconds to search for the peak correlation (default: 350). Only used when metric is "ccg".

Returns:

Pairwise similarity scores between

all unit pairs for each slice. Shape is (U, U, S).

lag_stack (PairwiseCompMatrixStack or None): Lag at which maximum

similarity occurs for each pair per slice. Shape is (U, U, S). None when metric is "sttc" (STTC has no lag).

av_corr (np.ndarray): Average similarity per slice across all unit

pairs in the lower triangle. Shape is (S,).

av_lag (np.ndarray or None): Average lag per slice. Shape is (S,).

None when metric is "sttc".

Return type:

corr_stack (PairwiseCompMatrixStack)

Notes

• Analogous to RateSliceStack.unit_to_unit_correlation but operates on raw spike trains instead of firing rate time series.

get_slice_to_slice_unit_comparison(metric='ccg', delt=20.0, bin_size=1.0, max_lag=350, min_spikes=2, min_frac=0.3, frac_active=None, n_jobs=-1)

Compute slice-to-slice similarity for each unit using spike-based metrics.

For each unit independently, compares its spike train across every pair of slices. Asks: “Does unit X fire in the same temporal pattern across repeated events?” Returns a PairwiseCompMatrixStack (S, S, U).

Parameters:

• metric (str) – Similarity metric to use. "ccg" for cross-correlogram on binned rasters (default), "sttc" for spike time tiling coefficient.

• delt (float) – STTC time window in milliseconds (default: 20.0). Only used when metric is "sttc".

• bin_size (float) – Bin size in milliseconds for the binary raster (default: 1.0). Only used when metric is "ccg".

• max_lag (float) – Maximum lag in milliseconds to search for the peak correlation (default: 350). Only used when metric is "ccg".

• min_spikes (int) – Minimum number of spikes in a slice for a unit to be considered valid in that slice (default: 2).

• min_frac (float) – Maximum fraction of slices that can be invalid before a unit’s average is set to NaN (default: 0.3).

• frac_active (np.ndarray or None) – Optional pre-computed fraction-active array of shape (U,) to override the internal per-unit validity check for computing averages. When provided, a unit’s average is set to NaN if frac_active[u] < (1 - min_frac). min_spikes still controls which individual slice pairs are computed. Compatible sources: SpikeSliceStack.compute_frac_active and SpikeData.get_frac_active (frac_per_unit output).

• n_jobs (int) – Number of threads for parallel computation. -1 uses all cores (default), 1 disables parallelism, None is serial.

Returns:

Pairwise similarity between all

slice pairs for each unit. Shape is (S, S, U).

all_lag (PairwiseCompMatrixStack or None): Lag at which maximum

similarity occurs for each slice pair per unit. Shape is (S, S, U). None when metric is "sttc".

av_corr (np.ndarray): Average similarity per unit across all valid

slice pairs. Shape is (U,).

av_lag (np.ndarray or None): Average lag per unit. Shape is (U,).

None when metric is "sttc".

Return type:

all_corr (PairwiseCompMatrixStack)

Notes

• Analogous to RateSliceStack.get_slice_to_slice_unit_corr_from_stack but operates on raw spike trains.

• Spike times within each slice are relative to the slice time origin (0-based or event-centered) for aligned comparison.

get_unit_timing_per_slice(timing='median', min_spikes=2)

Compute a representative spike time for each unit in each slice.

Returns a (U, S) matrix where entry [u, s] is the timing value (in milliseconds relative to the slice’s time origin) for unit u in slice s. For event-centered slices, t=0 is the event. Units with fewer than min_spikes spikes in a slice are marked NaN.

Parameters:

• timing (str) – Which spike time to extract per unit per slice. "median" (default), "mean", or "first" (onset latency).

• min_spikes (int) – Minimum number of spikes for a unit to count as active in a slice (default: 2).

Returns:

Array of shape (U, S) with timing

values in milliseconds relative to each slice’s time origin. NaN where the unit is inactive.

Return type:

timing_matrix (np.ndarray)

Notes

• Values are in milliseconds, not bin indices. This differs from RateSliceStack.get_unit_timing_per_slice which returns bin indices (suitable for direct indexing into the event stack). Both representations preserve rank order, so rank_order_correlation produces identical results either way.

• The returned matrix can be passed to rank_order_correlation to compute Spearman rank correlations between slice pairs, or used as input to order_units_across_slices for manual inspection of per-slice timing values.

rank_order_correlation(timing_matrix=None, timing='median', min_spikes=2, min_overlap=3, n_shuffles=100, min_overlap_frac=None, seed=1, n_jobs=-1)

Compute Spearman rank-order correlation of unit timing between all slice pairs.

For each pair of slices, only units active in both slices (non-NaN in both columns of the timing matrix) are included. If the overlap falls below the required minimum, the pair is set to NaN.

When n_shuffles > 0, the rank orders are shuffled n_shuffles times for each pair to build a null distribution, and the raw correlation is z-score normalised against it.

Parameters:

• timing_matrix (np.ndarray or None) – Array of shape (U, S) with timing values per unit per slice. NaN entries mark inactive units. Typically produced by get_unit_timing_per_slice. When None, computed automatically using timing and min_spikes.

• timing (str) – Which spike time to extract per unit per slice. "median" (default), "mean", or "first". Only used when timing_matrix is None.

• min_spikes (int) – Minimum spikes for activity (default: 2). Only used when timing_matrix is None.

• min_overlap (int) – Minimum number of units that must be active in both slices (default: 3).

• min_overlap_frac (float or None) – Minimum fraction of total units that must be active in both slices (default: None). When provided, the effective threshold is max(min_overlap, ceil(min_overlap_frac * U)).

• n_shuffles (int) – Number of shuffle iterations for z-scoring (default: 100). Set to 0 to return raw Spearman correlations. Values between 1 and 4 are rejected (minimum 5 required for a meaningful null distribution).

• seed (int or None) – Random seed for reproducibility of the shuffle (default: 1).

• n_jobs (int) – Number of threads for parallel computation. -1 uses all cores (default), 1 disables parallelism, None is serial.

Returns:

Spearman correlation matrix of

shape (S, S). When n_shuffles > 0, values are z-scores. When n_shuffles == 0, values are raw Spearman correlations.

av_corr (float): Average correlation (or z-score) across all valid

lower-triangle pairs.

overlap_matrix (PairwiseCompMatrix): Matrix of shape (S, S)

with fraction of units active in both slices.

Return type:

corr_matrix (PairwiseCompMatrix)

plot_aligned_slice_single_unit(unit_idx, ax=None, color_vals=None, color_label='', cmap='viridis', time_offset=0, xlabel='Rel. time (ms)', ylabel='Burst', x_range=None, vlines=None, show_colorbar=True, marker_size=20, font_size=None, style='scatter', invert_y=False, linewidths=0.5)

Plot a single unit’s spike times across all slices as a raster.

Extracts the spike train for unit_idx from every slice and delegates to plot_aligned_slice_single_unit().

Parameters:

• unit_idx (int) – Index of the unit to plot.

• ax (matplotlib.axes.Axes or None) – Target axes. If None, a new figure and axes are created.

• color_vals (np.ndarray or None) – Per-slice colour values.

• color_label (str) – Colorbar label.

• cmap (str) – Matplotlib colormap name.

• time_offset (float) – Value subtracted from every spike time before plotting. Slices from align_to_events are already event-centered (spike times in [-pre_ms, +post_ms]), so use the default time_offset=0. Only set a non-zero value when spike times are not already centered on the event.

• xlabel (str) – X-axis label.

• ylabel (str) – Y-axis label.

• x_range (tuple or None) – (xmin, xmax) for the x-axis.

• vlines (list[dict] or None) – Vertical reference lines. Each dict must contain 'x' and may optionally include 'color', 'linestyle', 'linewidth'.

• show_colorbar (bool) – Add a colorbar when color_vals is provided.

• marker_size (float) – Scatter marker size.

• font_size (int or None) – Font size for labels/ticks.

• style (str) – "scatter" for dot markers, "eventplot" for vertical line markers.

• invert_y (bool) – If True, first slice at top, last at bottom.

• linewidths (float) – Line width for eventplot markers.

Returns:

(fig, ax, sc) when ax is None, otherwise just sc.

sc is the scatter PathCollection (or None if no colour coding).

Return type:

result