I would like to find bursts among neuronal spikes.

20 views (last 30 days)
I have an array with neuronal interspike intervals (ISI). I would like to recognise bursts according to the following parameters: a burst starts if the ISI is lower than 6 ms, and ends if the ISI is above 9 ms. The recquirements of a burst is min 3 spikes, and 5 ms in length. There should be 20 ms between 2 bursts.
Can you help me which command should be used to find the bursting spikes? I need the number of bursts.

Accepted Answer

MarKf
MarKf on 5 Mar 2023
There are several features that you may want to take into account depending on your data and on what you want to do. The simplest thing would be to set up a threshold (say, 5 std above, ISI as a diff of spike indexes) that would not capture different types of spikes, tails, and other things, but still that might be enough for you. Have a look at some tools other people have developed for this, for example in fieldtrip (toolbox), here (this has minNspykes, duration and tails) and here (combination of multiple amplitude and spectral features, common approaches).
  4 Comments
MarKf
MarKf on 7 Mar 2023
Edited: MarKf on 7 Mar 2023
Actually I see now that the name ISIeg was misleading, sorry, it's not inter-spikes intervals. ISIeg was an example of indices of spikes (when they occur, if you had logical indices instead, i.e. zeros and a 1 every time there is a spike, you could just find them with find(logical_idx) but nevermind that now).
Your shared file is actual ISI probably, so what would become (unpadded) ISId in the previous script. Then it's the same:
ISIdat = load(websave('rd', "https://nl.mathworks.com/matlabcentral/answers/uploaded_files/1317005/ISI.mat"));
ISI = ISIdat.ISI; % this is just load('ISI.mat') file you shared above
spikes = cumsum(ISI); %what was ISIeg
bursts = zeros([sum(ISI), 1]); %or spikes(end)
bursts(spikes) = 0.5; % to plot bursts & interbursts
minI = 6; maxI = 9; minS = 3; lenB = 5; ibi = 20; %parameters
bibi = ibi*2; % to pad
ISId = [bibi ISI' bibi];
startb = ISId < minI;
contb = ISId < maxI;
starts = strfind([false, startb], [0 1]);
Nbursts = 0;
bbepv = -bibi; %initialise previous burst end
for sti = 1:numel(starts)
bst = starts(sti)-1;
bis = true; cou = 0;
while bis
cou = cou + 1; %next spike
bend = bst+cou; % cou(nter)+1 how many spikes
if ~contb(bend+1), bis = false; end % if ISIdiff > maxISI
end
bdur = spikes(bend)- spikes(bst); % burst duration
bint = spikes(bst)-bbepv; % since last burst
if cou+1 >= minS && bdur >= lenB && bint >= ibi %check "recquirements"
bbepv = spikes(bend); %in theory you don't have to recreate spikes indices
bursts(spikes(bst):bbepv) = 1; % you could do everything with ISI
Nbursts = Nbursts+1; %but the code was already here done, so I guess the above is left as an excercise?
end
end
plot(bursts), title(['N bursts = ' num2str(Nbursts)]),
ylim([0 1.5]), ylabel('interbursts = 0.5; bursts = 1')
Make sure the code does what it is supposed to do, with my quick 3 spikes example I created I knew where and how many spikes there were. Here there are many frequent interbusts, so maybe a smaller ibi would make more sense. Consider accepting the answer if that's what you were looking for or if it helped.
Bence Laczó
Bence Laczó on 12 Mar 2023
I tried your code but it found 10 bursts among the spikes. With manual counting I found only 4 bursts. I changed the recquirements in the 24th row from if cou+1 >= minS && bdur >= lenB && bint >= ibi to if cou >= minS && bdur >= lenB && bint >= ibi. Now it works fine. Does it make sense? Why did you write cou+1? I don't really understand it.

Sign in to comment.

More Answers (0)

Categories

Find more on Neural Simulation in Help Center and File Exchange

Products

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!