Finding peaks in EEG data

I know this question has been asked a number of times, but a thorough search did not help me find an answer and as I am very new to matlab I need some help. I need to count the number of peaks in a large body of EEG data. The code below allows me to plot what the data looks like (for 1 piece of data). I need to count the really long narrow spikes.

I know there are a bunch of peakfinder functions out there but I have no idea how to use them because they are all vector based and I don't know how to reference the vectors that are preprocessed in the script (fieldtrip). A simple threshold function wouldn't work because there is a lot of overall noise that is not related to what I want to count.

This is what the data looks like:

This is the script I have:

        % create the trial definition
        cfg=[];
        cfg.filename              = ['I:\EEG\',SubNum{subject},'\',FolderName, GroupName, SeqNum{subject}, '_', SubNum{subject}, '_Block', num2str(block), '.bdf'];
        cfg.dataset               = cfg.filename;

        cfg.continuous            = 'yes';

        data_org = ft_preprocessing(cfg)
        channel = 164;

        plot(data_org.time{1}, data_org.trial{1}(channel, :))
        xlabel('time (s)')
        ylabel('channel amplitude (uV)')
        legend(data_org.label(channel))

I don't know if it's asking too much, but I would really appreciate some guidance, I've been trying to use GUI's forever and it's no luck...but my coding skills are completely insufficient.

Answer 1

[Image Analyst]

You forgot to attach your data, so we can't try anything with it. Attach it if you want better help. Otherwise, if you have the Signal Processing Toolbox, you can try medfilt(), sgolayfilt(), and findpeaks(). It sounds like you tried findpeaks() and say it doesn't work for some reason, though you did not specify why it failed or what parameters you used for it.

For the Savitzky-Golay filter, use it to smooth the signal and remove spikes. This gives you like a background that you can then subtract from the original signal to give a "spikes-only" signal. Then you can use a global threshold. You can use the same basic process with the median filter, or even conv() for that matter. It looks like spikes that are more than about 20 are legitimate spikes and those less than 20 are just noise, but you can play with that number to set it optimally.

Answer 2

[Star Strider]

The first thing I would do is to detrend the baseline on a copy of your original data. That will make your task easier. [The polyfit (and polyval) functions are likely best for this.] You may have to do this more than once, since it’s best to use polynomial orders <10.

Digitally filtering your data might also be beneficial, at the risk of eliminating some of the data in the spikes. But if all you want to do is count them, an appropriately-designed bandpass filter might work for both baseline filtering and noise reduction, and would likely be better than simply detrending the data.

After preprocessing with baseline correction, bandpass filtering, or both, the Signal Processing Toolbox findpeaks function or a threshold function might be all that’s necessary.

AFTERTHOUGHT — Also consider using wavelets. They might be most appropriate for both denoising your signal and isolating the spikes of interest so you can count them.

[Image Analyst]

The Savitzky-Golay filter is like polyfit() except it does it just on data within a little window that slides along, so it will be better at "detrending" the data than using polyfit on the whole dataset to create just one big global quadratic or cubic. I have a demo if needed.

I would think that findpeaks() should be able to do a decent job without leveling the curve by subtracting the "background" smooth baseline signal if you give a good 'threshold' parameter. He said that he didn't know how to pass his signal into findpeaks(). I don't know why - seems pretty straightforward to me. Just pass in the 1D array that is the signal, and pass in a threshold. What's the big deal?