CofiringCoefficient

CofiringCoefficient - Compute the cofiring coefficient between two sets of intervals.

 Compute the cofiring coefficient of unit pairs between two sets of intervals
 (such as theta cycles and ripples). First, count the spikes emitted by each
 unit in each interval of the first set. Compute Pearson's correlation between
 the spike counts for each cell pair, and store the results in a vector. Then
 repeat this for the second set of intervals. The cofiring coefficient is
 defined as Pearson's correlation between the two resulting vectors. Thus,
 a high value indicates that cell pair correlations are similar between the two
 sets of intervals.

  USAGE

    [r,p] = CofiringCoefficient(spikes,intervals1,intervals2)

  INPUT

    spikes         two-column matrix of timestamps and unit IDs,
                   provided by <a href="matlab:help GetSpikeTimes">GetSpikeTimes</a> using 'output' = 'numbered'
    intervals1     two-column matrix of start and end times
    intervals2     two-column matrix of start and end times

  OUTPUT

    r              Pearson's r between the neuronal pair correlations in
                   intervals1 and in  intervals2 (cofiring coefficient).
    p              p-value for Pearson's test

0001 function [r,p] = CofiringCoefficient(spikes,intervals1,intervals2)
0002 
0003 %CofiringCoefficient - Compute the cofiring coefficient between two sets of intervals.
0004 %
0005 % Compute the cofiring coefficient of unit pairs between two sets of intervals
0006 % (such as theta cycles and ripples). First, count the spikes emitted by each
0007 % unit in each interval of the first set. Compute Pearson's correlation between
0008 % the spike counts for each cell pair, and store the results in a vector. Then
0009 % repeat this for the second set of intervals. The cofiring coefficient is
0010 % defined as Pearson's correlation between the two resulting vectors. Thus,
0011 % a high value indicates that cell pair correlations are similar between the two
0012 % sets of intervals.
0013 %
0014 %  USAGE
0015 %
0016 %    [r,p] = CofiringCoefficient(spikes,intervals1,intervals2)
0017 %
0018 %  INPUT
0019 %
0020 %    spikes         two-column matrix of timestamps and unit IDs,
0021 %                   provided by <a href="matlab:help GetSpikeTimes">GetSpikeTimes</a> using 'output' = 'numbered'
0022 %    intervals1     two-column matrix of start and end times
0023 %    intervals2     two-column matrix of start and end times
0024 %
0025 %  OUTPUT
0026 %
0027 %    r              Pearson's r between the neuronal pair correlations in
0028 %                   intervals1 and in  intervals2 (cofiring coefficient).
0029 %    p              p-value for Pearson's test
0030 %
0031 
0032 % Copyright (C) 2014-2018 by Ralitsa Todorova and Michaël Zugaro
0033 %
0034 % This program is free software; you can redistribute it and/or modify
0035 % it under the terms of the GNU General Public License as published by
0036 % the Free Software Foundation; either version 3 of the License, or
0037 % (at your option) any later version.
0038 
0039 if nargin < 3,
0040     error('Incorrect number of parameters (type ''help <a href="matlab:help CofiringCoefficient">CofiringCoefficient</a>'' for details).');
0041 end
0042 
0043 if isempty(spikes) || ~isdmatrix(spikes,'@2') || ~isivector(spikes(:,2)),
0044     error('Incorrect spikes (type ''help <a href="matlab:help CofiringCoefficient">CofiringCoefficient</a>'' for details).');
0045 end
0046 if isempty(intervals1) || ~isdmatrix(intervals1,'@2'),
0047     error('Incorrect intervals1 (type ''help <a href="matlab:help CofiringCoefficient">CofiringCoefficient</a>'' for details).');
0048 end
0049 if isempty(intervals2) || ~isdmatrix(intervals2,'@2'),
0050     error('Incorrect intervals2 (type ''help <a href="matlab:help CofiringCoefficient">CofiringCoefficient</a>'' for details).');
0051 end
0052 
0053 % Compute Pearson's r between pairs of units, resulting in a column vector, where each line is a pair of units.
0054 rPairs1 = reshape(CorrInIntervals(spikes,intervals1),[],1);
0055 rPairs2 = reshape(CorrInIntervals(spikes,intervals2),[],1);
0056 
0057 % Calculate Pearson's correlation between the pair correlations.
0058 % (If a unit did not fire within intervals1 or intervals2, all pairs with this unit must be ignored,
0059 % otherwise r would be NaN).
0060 ok = ~isnan(rPairs1) & ~isnan(rPairs2); 
0061 if sum(ok) == 0,
0062     r = NaN;
0063     return
0064 end
0065 if nargout > 1,
0066     [r,p] = corr(rPairs1(ok),rPairs2(ok));
0067 else
0068     r = corr(rPairs1(ok),rPairs2(ok));
0069 end
0070 
0071 end
0072 
0073 % ------------------------------- Helper function -------------------------------
0074 
0075 function C = CorrInIntervals(spikes,intervals);
0076 
0077 % Computes Pearson's correlations for the spikes of each pair of units within a set of intervals.
0078 % The spikes emitted by each unit are counted in each interval, then the two spike counts are correlated.
0079 % (Cii is set to NaN).
0080 
0081 n = max(spikes(:,2)); % number of units
0082 
0083 
0084 % List the interval number and unitID for each spike.
0085 [inIntervals,intervalIDs] = InIntervals(spikes(:,1),intervals);
0086 spikeList = [intervalIDs(inIntervals) spikes(inIntervals,2)];
0087 
0088 % Generate the spike count matrix for each unit (column) and each interval (row).
0089 spikeCountMatrix = Accumulate(spikeList,1,[n size(intervals,1)]);
0090 
0091 % Compute the correlation matrix C. Because C is symmetrical (Cij=Cji), keep only upper triangular elements
0092 % in order to remove duplicates, so that each pair is counted only once.
0093 C = corr(spikeCountMatrix);
0094 C(logical(tril(ones(size(C))))) = NaN;
0095 
0096 end