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: AnEvent
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: ANeuralData
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.
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)
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 anEvent
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.
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
andbrowseSequences
should be modified/combined. The newbrowseEvents
should have signaturebrowseEvents(cellEvents, <nd>)
and only plot the spike-raster plot if noNeuralData
object is provided. The newbrowseSequences
should bebrowseSequences(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.