[Eeglablist] PSD on scored PSG data

[Makoto Miyakoshi]

Dear Mohith,

> Thank you for sharing the code but perhaps I am not making myself clear
enough that for the PSG data that I acquired and which have been event
coded into different stages of sleep I am required to run a frequency
domain analysis of visualizing as well as obtaining the mean power spectrum
for specific frequency bands like alpha, delta, theta etc FOR EACH SLEEP
STAGE. So running a spectrogram may not be the appropriate way to visualize
how channel or component activity between participants differ in terms of
power spectra within each sleep stage (e.g. REM vs. NREM 1 vs. NREM2 VS.
NREM3).

The obtained spectrogram can be also segmented into blocks, if you have
annotations, and you can take a block average.

> According to what I read upon this topic, MATLAB has functions like
pwelch or fft to run such sort of analysis and in EEGLAB, I think its
pop_spectopo instead of spectrogram (since I don't need info about the time
domain at least for this part of the analysis), please correct me if I am
wrong.

In the example below, you can see the PSD time in 'times' and freqs in
'freqs'.
 [~, freqs, times, PSD]      = spectrogram(EEG.icaact(icIdx,:), EEG.srate,
round(EEG.srate/2), freqBins, EEG.srate);
If you know your Stage 1 is from 568.9 s to 668.2 sec,
startIdx = find(times>568.9,1);
endIdx = find(times>668.2,1)-1;
stage1MeanPower = squeeze(mean(PSD(:, startIdx:endIdx, icIdx), 2));
to obtain the Stage1 PSD (of a given independent component by icIdx).

Makoto



On Sat, Jan 12, 2019 at 7:21 PM Makoto Miyakoshi <mmiyakoshi at ucsd.edu>
wrote:

> Dear Mohith,
>
> Oops, I pasted a wrong one. Here is the right version, sorry for
> confusion. Please discard the last one and use this one. Note the second
> input arguments to spectrogram() is EEG.srate which means 1-s sliding
> window is being used.
>
> Makoto
>
>
>
> % Define the freq bins.
> freqBins = 1:0.01:40;
>
> % Create an empty PSD tensor (time & freq & IC).
> PSD_tensor = [];
>
> % Compute the PSD tensor.
> for icIdx = 1:size(EEG.icaact,1)
>
>     % Compute short-term Fourier Transform for 1-s window, 1-50 Hz, 0.1 Hz
> bin.
>     if icIdx == 1;
>         [~, freqs, times, firstPSD]      =
> spectrogram(EEG.icaact(icIdx,:), EEG.srate, round(EEG.srate/2), freqBins,
> EEG.srate);
>         PSD_tensor = zeros(length(freqs), size(firstPSD,2),
> size(EEG.icaact,1));
>         PSD_tensor(:,:,icIdx) = firstPSD;
>     else
>         [~, ~, ~, PSD_tensor(:,:,icIdx)] =
> spectrogram(EEG.icaact(icIdx,:), EEG.srate, round(EEG.srate/2), freqBins,
> EEG.srate);
>     end
> end
>
> On Thu, Jan 10, 2019 at 6:53 PM Makoto Miyakoshi <mmiyakoshi at ucsd.edu>
> wrote:
>
>> Dear Mohith,
>>
>> See my latest response to Panos (pasted below).
>> In EEGLAB, your time-series data are stored under EEG.data as channels x
>> time. So, if you want to use spectrogram(), it will be something like
>> this...see below. This is from my recent work so it works. Note I'm working
>> on EEG.icaact which is one of ICA's results. Replace it with EEG.data if
>> you want to work on sensor level.
>>
>> Makoto
>>
>>
>>
>> % Define the freq bins.
>> freqBins = 1:0.01:50;
>>
>> % Create an empty PSD tensor (time & freq & IC).
>> PSD_tensor = [];
>>
>> % Compute the PSD tensor.
>> for icIdx = 1:size(EEG.icaact,1)
>>
>>     % Compute short-term Fourier Transform for 1-s window, 1-50 Hz, 0.1
>> Hz bin.
>>     if icIdx == 1;
>>         [~, freqs, times, firstPSD]      =
>> spectrogram(EEG.icaact(icIdx,:), ceil(EEG.xmax), [], freqBins, EEG.srate);
>>         PSD_tensor = zeros(length(freqs), size(firstPSD,2),
>> size(EEG.icaact,1));
>>         PSD_tensor(:,:,icIdx) = firstPSD;
>>     else
>>         [~, ~, ~, PSD_tensor(:,:,icIdx)] =
>> spectrogram(EEG.icaact(icIdx,:), ceil(EEG.xmax), [], freqBins, EEG.srate);
>>     end
>> end
>>
>> Makoto
>>
>>
>>
>> On Thu, Jan 10, 2019 at 3:57 PM Makoto Miyakoshi <mmiyakoshi at ucsd.edu>
>> wrote:
>>
>>> Dear Panos,
>>>
>>> > 1) In addition to calculating the average absolute power (as your
>>> script nicely shows), I was also interested in calculating the average
>>> absolute (and relative) power at binned time intervals (e.g. avg power
>>> between 0-1sec, avg power between 1-2sec, etc) within the dataset. I tried
>>> to use the "spectra" output from spectopo but from what I gather it comes
>>> up with [(sampling rate)/2 + 1] points rather that one power-spectral point
>>> per timepoint. How would you recommend that I proceed?
>>>
>>> You can repeatedly apply EEGLAB spectopo() function to perform hand-made
>>> short-term Fourier transform (STFT), but alternatively you might want to
>>> use either EEGLAB newtimef() or Matlab spectrogram() function (the latter
>>> may require some additional Toolbox). The output will be frequency x time
>>> matrix. The interval of time bins needs to be calcualted. Basically,
>>> {(length of data) - (sliding window length)}/(number of steps) gives you
>>> the interval (step size). Adjust the (number of steps) so that you can
>>> obtain the desired interval.
>>>
>>> > 2) Is there a way to display how do topographic maps (scalp heat maps)
>>> change with time (I'm able to see how they change with different
>>> frequencies but I was interested in seeing how they also change with time)?
>>> Would the function timtopo be the best way to do that?
>>>
>>> See this wiki page.
>>>
>>> https://sccn.ucsd.edu/wiki/Chapter_02:_Writing_EEGLAB_Scripts#Creating_a_scalp_map_animation
>>>
>>> > 3) A more general question:  If I write a matlab script that I would
>>> like to apply on a bunch of datasets (which in my case are just epochs of
>>> different lengths that I have extracted from my original dataset), should I
>>> put all said datasets  (which I have already pre-processed and applied ICA
>>> on) in a STUDY set and then apply the script there, or should I just write
>>> a for loop in matlab and apply the script in each individual dataset? In
>>> other words, does the STUDY set offer an advantage in this case?  (I
>>> apologize for the potential triviality of this one!)
>>>
>>> If you are a beginner, it is always a good idea to make things as simple
>>> as possible. I recommend you organize your own code to loop the
>>> single-subject process for all the subjects. After all, that's the only to
>>> learn the process!
>>>
>>> Makoto
>>>
>>> On Wed, Jan 9, 2019 at 11:47 AM Fotiadis, Panagiotis <
>>> Panagiotis.Fotiadis at pennmedicine.upenn.edu> wrote:
>>>
>>>> Hi Makoto,
>>>>
>>>>
>>>> Thank you for the really great advice! The two links you provided are
>>>> extremely helpful.
>>>>
>>>>
>>>> I had a few follow-up questions:
>>>>
>>>> 1) In addition to calculating the average absolute power (as your
>>>> script nicely shows), I was also interested in calculating the average
>>>> absolute (and relative) power at binned time intervals (e.g. avg power
>>>> between 0-1sec, avg power between 1-2sec, etc) within the dataset. I tried
>>>> to use the "spectra" output from spectopo but from what I gather it comes
>>>> up with [(sampling rate)/2 + 1] points rather that one power-spectral point
>>>> per timepoint. How would you recommend that I proceed?
>>>>
>>>>
>>>> 2) Is there a way to display how do topographic maps (scalp heat maps)
>>>> change with time (I'm able to see how they change with different
>>>> frequencies but I was interested in seeing how they also change with time)?
>>>> Would the function timtopo be the best way to do that?
>>>>
>>>>
>>>> 3) A more general question:  If I write a matlab script that I would
>>>> like to apply on a bunch of datasets (which in my case are just epochs of
>>>> different lengths that I have extracted from my original dataset), should I
>>>> put all said datasets  (which I have already pre-processed and applied ICA
>>>> on) in a STUDY set and then apply the script there, or should I just write
>>>> a for loop in matlab and apply the script in each individual dataset? In
>>>> other words, does the STUDY set offer an advantage in this case?  (I
>>>> apologize for the potential triviality of this one!)
>>>>
>>>>
>>>> Thank you again in advance for your time and help!
>>>>
>>>>
>>>> Best,
>>>>
>>>> Panos
>>>>
>>>>
>>>> Panagiotis Fotiadis
>>>>
>>>> PhD Student | Neuroscience Graduate Group
>>>>
>>>> Perelman School of Medicine, University of Pennsylvania
>>>> ------------------------------
>>>> *From:* Makoto Miyakoshi <mmiyakoshi at ucsd.edu>
>>>> *Sent:* Monday, January 7, 2019 2:48:37 PM
>>>> *To:* Fotiadis, Panagiotis
>>>> *Cc:* eeglablist at sccn.ucsd.edu
>>>> *Subject:* [External] Re: [Eeglablist] Frequency-time spectrogram
>>>> deconstruction
>>>>
>>>> Dear Panos,
>>>>
>>>> Welcome to the time-frequency world.
>>>>
>>>> > Would I just need to bandpass filter my post-processed EEG signal to
>>>> each frequency range of interest (i.e., alpha: 8-12Hz etc) and then plot
>>>> the remaining EEG signal over time, or is there another way to do this?
>>>>
>>>> That's one way to go. Nothing is wrong with that!
>>>>
>>>> More convenient and established way to go is to perform time-frequency
>>>> transform using short-term Fourier transform or Wavelet transform. Google
>>>> EEGLAB time-frequency and you'll find many of our past workshop materials.
>>>> For example, see Slide 21 of this file
>>>>
>>>> https://sccn.ucsd.edu/mediawiki/images/a/a6/C2_A3_Time-frequencyDecAndAdvancedICAPracticum_updateJan2017.pdf
>>>>
>>>> You can also obtain bin-mean values from power spectral density. See
>>>> below.
>>>>
>>>> https://sccn.ucsd.edu/wiki/Makoto's_useful_EEGLAB_code#How_to_extract_EEG_power_of_frequency_bands
>>>>
>>>> Makoto
>>>>
>>>> On Mon, Jan 7, 2019 at 1:34 AM Fotiadis, Panagiotis <
>>>> Panagiotis.Fotiadis at pennmedicine.upenn.edu> wrote:
>>>>
>>>> Hello,
>>>>
>>>>
>>>> I am fairly new to EEGLab and I had a question concerning the
>>>> deconstruction of my EEG signal into its alpha/beta/theta/delta
>>>> sub-components:
>>>>
>>>>
>>>> After pre-processing some subjects with EEG data from 128 channels and
>>>> performing ICA (using runica), I used eeglab and chronux to plot the
>>>> power/frequency and frequency/time spectrograms of several epochs of
>>>> interest.
>>>>
>>>>
>>>> Is there a way to extract the alpha/beta/theta/delta frequencies of
>>>> those epochs and quantify when they occur in time? I can visualize when
>>>> each type of neuronal oscillation occurs by looking at the overall
>>>> frequency/time spectrogram, but I was wondering whether there was a more
>>>> robust way to actually plot each type of oscillation separately and/or
>>>> quantify when it occurs.
>>>>
>>>>
>>>> Would I just need to bandpass filter my post-processed EEG signal to
>>>> each frequency range of interest (i.e., alpha: 8-12Hz etc) and then plot
>>>> the remaining EEG signal over time, or is there another way to do this?
>>>>
>>>>
>>>> Thank you in advance!
>>>>
>>>>
>>>> Best,
>>>>
>>>> Panos
>>>> _______________________________________________
>>>> Eeglablist page: http://sccn.ucsd.edu/eeglab/eeglabmail.html
>>>> To unsubscribe, send an empty email to
>>>> eeglablist-unsubscribe at sccn.ucsd.edu
>>>> For digest mode, send an email with the subject "set digest mime" to
>>>> eeglablist-request at sccn.ucsd.edu
>>>>
>>>>
>>>>
>>>> --
>>>> Makoto Miyakoshi
>>>> Swartz Center for Computational Neuroscience
>>>> Institute for Neural Computation, University of California San Diego
>>>>
>>>
>>>
>>> --
>>> Makoto Miyakoshi
>>> Swartz Center for Computational Neuroscience
>>> Institute for Neural Computation, University of California San Diego
>>>
>>
>>
>> --
>> Makoto Miyakoshi
>> Swartz Center for Computational Neuroscience
>> Institute for Neural Computation, University of California San Diego
>>
>> On Thu, Jan 10, 2019 at 6:34 PM VARMA, MOHITH MUKUND <
>> mohith96 at connect.hku.hk> wrote:
>>
>>> Dear all,
>>>
>>> I am trying to conduct PSD (Power Spectrum Density) analysis for each of
>>> the sleep stages across different frequency bands. I have scored the PSG
>>> data (edf format) on a Python toolbox called Sleep and obtained the
>>> hypnogram that I can use as time codes for each sleep stage. I imported
>>> this time code file as event file (it contains info on the latency and type
>>> of sleep stage) and now the thing I am stuck on is how to run PSD (possibly
>>> using spectopo) for the same event type across different frequency band. I
>>> think I cannot use the GUI to run this part so your help for setting up the
>>> parameters in the MATLAB command would be really helpful! My data's
>>> sampling rate is 250 Hz and FFT window size is 4 sec and other paramters I
>>> can go with the default settings provided by EEGLAB. I hope my question is
>>> clear enough, please let me know if you need further information.
>>>
>>> Regards,
>>>
>>> --
>>> Mohith M. Varma (Mo)
>>> Graduate Research Assistant
>>>
>>> Social & Cognitive Neuroscience Laboratory
>>> Department of Psychology
>>> Faculty of Social Sciences
>>> The University of Hong Kong
>>> Tel: (+852) 52622875
>>> Email: mohith96 at connect.hku.hk
>>> _______________________________________________
>>> Eeglablist page: http://sccn.ucsd.edu/eeglab/eeglabmail.html
>>> To unsubscribe, send an empty email to
>>> eeglablist-unsubscribe at sccn.ucsd.edu
>>> For digest mode, send an email with the subject "set digest mime" to
>>> eeglablist-request at sccn.ucsd.edu
>>
>>
>>
>> --
>> Makoto Miyakoshi
>> Swartz Center for Computational Neuroscience
>> Institute for Neural Computation, University of California San Diego
>>
>
>
> --
> Makoto Miyakoshi
> Swartz Center for Computational Neuroscience
> Institute for Neural Computation, University of California San Diego
>


-- 
Makoto Miyakoshi
Swartz Center for Computational Neuroscience
Institute for Neural Computation, University of California San Diego
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://sccn.ucsd.edu/pipermail/eeglablist/attachments/20190114/a5a9c18f/attachment.html>