I wish to compute the “delta power” of EEG recorded over the course of a night’s sleep. From the raw EEG signal (in microvolts), I’m attempting to do the short-time Fourier transform (STFT) in small windows of the raw signal, then analyze the outputs in the range of 0.5-4Hz (delta).
I have lots of questions about 1) how STFTs should be ideally computed; 2) how Matlab computes STFTs; and 3) how to average results across frequency bands and time ranges of interest.
I apologize for the length of this post. It is clear that I should take a course in signal processing, but life is short, I’m in the middle of my career, and there is no time to learn things from 1st principles. Any help relieving my ignorance about all this will be greatly appreciated.
As background, I am running Matlab 126.96.36.1992 with the Signal Processing Toolbox on a WinXP machine (64 bit).
1) My understanding is that I could take a spectrogram of my raw EEG signal as a way to implement STFTs. Does this seem reasonable? As opposed to snipping out small time windows of the raw signal, hanning windowing them, computing FFT, calculating power of output, moving on to next window, lather rinse repeat, etc…
2) For some reason, I have two spectrogram functions: “spectrogram” and “specgram”. The latter may be 3rd party. When I try to call the spectrogram, it gives me the error: Undefined function or method 'gencoswin' in ‘hamming’ at 14. Is one of specgram / spectrogram not included in the signal processing toolbox? Is this a licensing error? Is there any great difference between these two functions?
3) The few papers I’ve read on the process of extracting delta power breaks the signal into 4s windows (=nfft), computes the STFT, then averages across several of those windows (e.g., to get avg power in 20s windows). Some authors also average across frequency bin in order to get a single estimate of power in the range 0.5-4Hz. It is usually not specified how this averaging is done: since abs, square and log are non-linear, it matters when I average across frequency and time bins. Any suggestions at what step averaging should be done (e.g., should I average the log PSD or take the log PSD of the average?)
4) Related to 3), is the output b from specgram (or s in spectrogram) in squared units?
5) Right now, I’m using a window of the size of nfft. In this case, I assume there can be no overlap, is that true? Is there a problem not overlapping the signal?
6) Let's say I take a 20s snippet of my raw signal. Then I calculate a spectrogram with a window size (nfft) of 4s. I would expect there to be 5 time bins and Nyquist*4 (+1) frequency bins. I do get this result (e.g., f = Nyquist*4+1; t = 5), but when I plot the output (as surf), it goes from 0-16 sec, not 0-20 sec. That is, it is plotting 4 time windows from 0-4,4-8,8-12,12-16s. If I increase my time resolution (but worsen my spectral resolution), I get something closer to 20s windows. Is there an ideal trade-off? If I go with nfft=4s, then should I start the next 20s time window one 4s window behind, in order to make up for the lost window from 16-20s?
7) In older papers, the window is a boxcar the length of nfft. I think specgram takes a hanning window. Is there a way to specify boxcar for specgram?
8) Finally, another approach would be to do all this in the time domain: generate a, say, Chebyshev Type I bandpass filter between 0.5 and 4Hz; filter the raw EEG signal with these coeffs (using filtfilt to avoid phase shift errors); then simply compute 1og power on the filtered values of the form 10*log10(x2.^2 + 1) (the +1 is to remove fractions in the signal, which will produce negative values after logging). Should this produce the same thing as the spectrogram version above?
Thanks for the help