This code is designed to provide tools for working with neuronal sequences. The project started originally as a set of functions to detect sharp-wave ripples in the hippocampus. The code breaks down into three basic pieces:

  • The sequences folder: At the core of this project, the aim was to analyze neuronal sequences. Our concept of a sequence is simply a list of spikes. Because we discard timing information from the spike trains, the techniques implemented here are agnostic of the timescale on which such a sequence was actually recorded.
  • The Event class: An Event is a convenience class that was designed (1) to store a sequence and additional information about its origin in a recording and (2) to bundle together functions for working with sequences.
  • The NeuralData class: A NeuralData object is our interface to the raw data from a neural recording. This class provides functions for retrieving recording information, extracting LFPs, and detecting various event types, among some other convenience functions.

The sections below describe each of these pieces in more detail along with some basic examples of how to use this project. This project is dependent on some additional functionality found in the matlab-incremented project.

Sequences

Since a sequence as we have it defined is just a list of spikes (i.e., neurons), each sequence represented in this project as a vector. As such, many MATLAB functions for vectors can come in handy. For instance, the collection of neurons that are active in a sequence s can be found with unique(s).

s = [1, 2, 1, 3, 2, 2, 4, 3, 3, 4]
unique(s)

A few sequence-specific functions have been written where there was no equivalent for vectors. The spike-count vector can be computed with spikecount, and the index sets where neurons spike can be extracted with spikesets.

spikecount(s)
spikesets(s)

To visualize s, create a spike-raster plot.

spikeraster(s)

The center-of-mass ordering of a sequence s can be computed with sortneurons(s).

sortneurons(s)

The bias-count and skew-bias matrices of s can be computed with biascount(s) and skewbias(s), respectively. Both of these functions return a sparse matrix.

full(biascount(s))
full(skewbias(s))

The spikes of a sequence can be permuted to be in a random order using shuffle.

u = shuffle(s)

To compute the correlation between sequences, use correlation.

correlation(s, u)

For a list L1 of sequences (i.e., a vector cell array), use activity(L1) to return a matrix in which row i is an indicator vector for the set of neurons that are active in L{i}.

randseq = @(~) randi(5, 1, 3);
L1 = arrayfcn(randseq, zeros(3, 1));
L1{:}
activity(L1)

With a second list L2 of sequences, the number of neurons that each pair of sequences has in common can be found with coactivity(L1, L2).

L2 = arrayfcn(randseq, zeros(3, 1));
L2{:}
coactivity(L1, L2)

Events

When dealing with sequences in real data, it is often necessary for purposes of analysis to know more about each sequence. In concept, an event is all of the information about a recording that occurs within some time window; as such, the Event object class stores this time window. However, extracting information from a NeuralData object can be time-consuming; so this class also stores some additional information. The members of this class are listed below:

window
This parameter stores the time window during which an event occurs. It is stored as a vector of the form [start, end]. For consistency regardless of the sampling rate of a recording, this is stored with units of seconds.
spikes
This parameter stores the list of spikes that occur during an event. That is, spikes is the neuronal sequence contained within an Event object.
times
This parameter stores the times corresponding to the event’s spikes.
type
This parameter allows the user to label each event as some specific type. A type can be an arbitrary string.

Though an Event can be created manually, it is often most useful in that several methods of the NeuralData class return (lists of) events.

Interfacing with real data

The code in the NeuralData class was originally designed to provide functionality for detecting sharp-wave ripples in LFP signals.

  • TODO: The current setup requires a .dat file with a specific naming convention. It also requires a *_BehavElecDataLFP.mat file and some of its kin. This should be updated so that these extra files are not necessary. In the meantime, mention these sisues and direct the user to the github page. As an intermediate step, perhaps add a function for generating these files. Make the .dat file optional if a cache directory exists. Use .res (etc.) files instead of the *_BehavElecDataLFP.mat file.

For the sake of example, suppose that the directory A111-20150701 contains the files A111-20150701.dat and A111-20150701_BehavElecDataLFP.mat. Also, assume that cache is a directory. We can create a NeuralData object by specifying the data and cache directories.

nd = NeuralData('A111-20150701', 'cache');

We can now extract various kinds of information from the object nd. To detect SWR events, three channels need to be selected by hand using an external tool. Once the channels have been selected, load them with loadChannels.

loadChannels(nd)

The SWR events can now be detected with detectRipples.

cellRipples = detectRipples(nd);

The output of this method is a cell array of Event objects. These events can be visualized with browseEvents.

browseEvents(nd, cellRipples)
  • TODO: The functions browseEvents and browseSequences should be modified/combined. The new browseEvents should have signature browseEvents(cellEvents, <nd>) and only plot the spike-raster plot if no NeuralData object is provided. The new browseSequences should be browseSequences(cellSeqs) = browseEvents(cellfcn(@seq2evt, cellSeq)).

The locations of the animal in the track at a particular list of times can be retrieved with the methods getLocationsAtTimes. For instance, we might want to know where ripples occurred in the track. To do this, we will need a single time to represent each event. Let’s choose the starting time of the events.

vTimes = cellfun(@startTime, cellRipples);
mtxLocs = getLocationsAtTimes(nd, vTimes);

We can now visualize where the ripples occur on the track.

vWinX = minmax(mtxLocs(:, 1));
vWinY = minmax(mtxLocs(:, 2));
dStd = 0.1 * min(diff(vWinX), diff(vWinY));
nResolution = 1000;

imagesc(psth2d(mtxLocs, vWinX, vWinY, dStd, nResolution))

Many other NeuralData methods exist for extraction of recording information and event information.