Compute source power spectral density (PSD) in a label#

Returns an STC file containing the PSD (in dB) of each of the sources within a label.

# Authors: Alexandre Gramfort <alexandre.gramfort@inria.fr>
#
# License: BSD-3-Clause
# Copyright the MNE-Python contributors.
import matplotlib.pyplot as plt

import mne
from mne import io
from mne.datasets import sample
from mne.minimum_norm import compute_source_psd, read_inverse_operator

print(__doc__)

Set parameters

data_path = sample.data_path()
meg_path = data_path / "MEG" / "sample"
raw_fname = meg_path / "sample_audvis_raw.fif"
fname_inv = meg_path / "sample_audvis-meg-oct-6-meg-inv.fif"
fname_label = meg_path / "labels" / "Aud-lh.label"

# Setup for reading the raw data
raw = io.read_raw_fif(raw_fname, verbose=False)
events = mne.find_events(raw, stim_channel="STI 014")
inverse_operator = read_inverse_operator(fname_inv)
raw.info["bads"] = ["MEG 2443", "EEG 053"]

# picks MEG gradiometers
picks = mne.pick_types(
    raw.info, meg=True, eeg=False, eog=True, stim=False, exclude="bads"
)

tmin, tmax = 0, 120  # use the first 120s of data
fmin, fmax = 4, 100  # look at frequencies between 4 and 100Hz
n_fft = 2048  # the FFT size (n_fft). Ideally a power of 2
label = mne.read_label(fname_label)

stc = compute_source_psd(
    raw,
    inverse_operator,
    lambda2=1.0 / 9.0,
    method="dSPM",
    tmin=tmin,
    tmax=tmax,
    fmin=fmin,
    fmax=fmax,
    pick_ori="normal",
    n_fft=n_fft,
    label=label,
    dB=True,
)

stc.save("psd_dSPM", overwrite=True)
320 events found on stim channel STI 014
Event IDs: [ 1  2  3  4  5 32]
Reading inverse operator decomposition from /home/circleci/mne_data/MNE-sample-data/MEG/sample/sample_audvis-meg-oct-6-meg-inv.fif...
    Reading inverse operator info...
    [done]
    Reading inverse operator decomposition...
    [done]
    305 x 305 full covariance (kind = 1) found.
    Read a total of 4 projection items:
        PCA-v1 (1 x 102) active
        PCA-v2 (1 x 102) active
        PCA-v3 (1 x 102) active
        Average EEG reference (1 x 60) active
    Noise covariance matrix read.
    22494 x 22494 diagonal covariance (kind = 2) found.
    Source covariance matrix read.
    22494 x 22494 diagonal covariance (kind = 6) found.
    Orientation priors read.
    22494 x 22494 diagonal covariance (kind = 5) found.
    Depth priors read.
    Did not find the desired covariance matrix (kind = 3)
    Reading a source space...
    Computing patch statistics...
    Patch information added...
    Distance information added...
    [done]
    Reading a source space...
    Computing patch statistics...
    Patch information added...
    Distance information added...
    [done]
    2 source spaces read
    Read a total of 4 projection items:
        PCA-v1 (1 x 102) active
        PCA-v2 (1 x 102) active
        PCA-v3 (1 x 102) active
        Average EEG reference (1 x 60) active
    Source spaces transformed to the inverse solution coordinate frame
Not setting metadata
70 matching events found
No baseline correction applied
Created an SSP operator (subspace dimension = 3)
3 projection items activated
Considering frequencies 4 ... 100 Hz
Preparing the inverse operator for use...
    Scaled noise and source covariance from nave = 1 to nave = 1
    Created the regularized inverter
    Created an SSP operator (subspace dimension = 3)
    Created the whitener using a noise covariance matrix with rank 302 (3 small eigenvalues omitted)
    Computing noise-normalization factors (dSPM)...
[done]
Picked 305 channels from the data
Computing inverse...
    Eigenleads need to be weighted ...
Reducing data rank 33 -> 33
Using hann windowing on at most 70 epochs

  0%|          |  : 0/70 [00:00<?,       ?it/s]
  1%|▏         |  : 1/70 [00:00<00:06,    9.98it/s]
  3%|▎         |  : 2/70 [00:00<00:04,   14.72it/s]
  4%|▍         |  : 3/70 [00:00<00:03,   19.20it/s]
  6%|▌         |  : 4/70 [00:00<00:02,   22.36it/s]
  7%|▋         |  : 5/70 [00:00<00:02,   24.87it/s]
  9%|▊         |  : 6/70 [00:00<00:02,   25.99it/s]
 10%|█         |  : 7/70 [00:00<00:02,   27.79it/s]
 11%|█▏        |  : 8/70 [00:00<00:02,   29.32it/s]
 13%|█▎        |  : 9/70 [00:00<00:01,   30.78it/s]
 14%|█▍        |  : 10/70 [00:00<00:01,   32.11it/s]
 16%|█▌        |  : 11/70 [00:00<00:01,   33.14it/s]
 17%|█▋        |  : 12/70 [00:00<00:01,   32.74it/s]
 19%|█▊        |  : 13/70 [00:00<00:01,   33.76it/s]
 20%|██        |  : 14/70 [00:00<00:01,   34.70it/s]
 21%|██▏       |  : 15/70 [00:00<00:01,   35.24it/s]
 23%|██▎       |  : 16/70 [00:00<00:01,   35.88it/s]
 24%|██▍       |  : 17/70 [00:00<00:01,   36.37it/s]
 26%|██▌       |  : 18/70 [00:00<00:01,   36.27it/s]
 27%|██▋       |  : 19/70 [00:00<00:01,   36.85it/s]
 29%|██▊       |  : 20/70 [00:00<00:01,   37.29it/s]
 30%|███       |  : 21/70 [00:00<00:01,   37.88it/s]
 31%|███▏      |  : 22/70 [00:00<00:01,   38.37it/s]
 33%|███▎      |  : 23/70 [00:00<00:01,   38.07it/s]
 34%|███▍      |  : 24/70 [00:00<00:01,   38.53it/s]
 36%|███▌      |  : 25/70 [00:00<00:01,   38.99it/s]
 37%|███▋      |  : 26/70 [00:00<00:01,   39.15it/s]
 39%|███▊      |  : 27/70 [00:00<00:01,   39.62it/s]
 40%|████      |  : 28/70 [00:00<00:01,   39.90it/s]
 41%|████▏     |  : 29/70 [00:00<00:01,   39.30it/s]
 43%|████▎     |  : 30/70 [00:00<00:01,   39.66it/s]
 44%|████▍     |  : 31/70 [00:00<00:00,   40.04it/s]
 46%|████▌     |  : 32/70 [00:00<00:00,   40.40it/s]
 47%|████▋     |  : 33/70 [00:00<00:00,   40.64it/s]
 49%|████▊     |  : 34/70 [00:00<00:00,   40.45it/s]
 50%|█████     |  : 35/70 [00:00<00:00,   40.71it/s]
 51%|█████▏    |  : 36/70 [00:00<00:00,   41.02it/s]
 53%|█████▎    |  : 37/70 [00:00<00:00,   41.24it/s]
 54%|█████▍    |  : 38/70 [00:00<00:00,   41.26it/s]
 56%|█████▌    |  : 39/70 [00:01<00:00,   41.38it/s]
 57%|█████▋    |  : 40/70 [00:01<00:00,   40.72it/s]
 59%|█████▊    |  : 41/70 [00:01<00:00,   40.97it/s]
 60%|██████    |  : 42/70 [00:01<00:00,   41.26it/s]
 61%|██████▏   |  : 43/70 [00:01<00:00,   41.32it/s]
 63%|██████▎   |  : 44/70 [00:01<00:00,   41.54it/s]
 64%|██████▍   |  : 45/70 [00:01<00:00,   41.60it/s]
 66%|██████▌   |  : 46/70 [00:01<00:00,   41.25it/s]
 67%|██████▋   |  : 47/70 [00:01<00:00,   41.30it/s]
 69%|██████▊   |  : 48/70 [00:01<00:00,   41.55it/s]
 70%|███████   |  : 49/70 [00:01<00:00,   41.82it/s]
 71%|███████▏  |  : 50/70 [00:01<00:00,   41.98it/s]
 73%|███████▎  |  : 51/70 [00:01<00:00,   41.59it/s]
 74%|███████▍  |  : 52/70 [00:01<00:00,   41.80it/s]
 76%|███████▌  |  : 53/70 [00:01<00:00,   41.97it/s]
 77%|███████▋  |  : 54/70 [00:01<00:00,   42.10it/s]
 79%|███████▊  |  : 55/70 [00:01<00:00,   42.31it/s]
 80%|████████  |  : 56/70 [00:01<00:00,   42.47it/s]
 81%|████████▏ |  : 57/70 [00:01<00:00,   41.92it/s]
 83%|████████▎ |  : 58/70 [00:01<00:00,   42.20it/s]
 84%|████████▍ |  : 59/70 [00:01<00:00,   42.45it/s]
 86%|████████▌ |  : 60/70 [00:01<00:00,   42.62it/s]
 87%|████████▋ |  : 61/70 [00:01<00:00,   42.74it/s]
 89%|████████▊ |  : 62/70 [00:01<00:00,   42.79it/s]
 90%|█████████ |  : 63/70 [00:01<00:00,   42.47it/s]
 91%|█████████▏|  : 64/70 [00:01<00:00,   42.65it/s]
 93%|█████████▎|  : 65/70 [00:01<00:00,   42.93it/s]
 94%|█████████▍|  : 66/70 [00:01<00:00,   42.98it/s]
 96%|█████████▌|  : 67/70 [00:01<00:00,   43.13it/s]
 97%|█████████▋|  : 68/70 [00:01<00:00,   42.75it/s]
 99%|█████████▊|  : 69/70 [00:01<00:00,   42.94it/s]
100%|██████████|  : 70/70 [00:01<00:00,   43.11it/s]
100%|██████████|  : 70/70 [00:01<00:00,   40.37it/s]
Writing STC to disk...
[done]

View PSD of sources in label

plt.plot(stc.times, stc.data.T)
plt.xlabel("Frequency (Hz)")
plt.ylabel("PSD (dB)")
plt.title("Source Power Spectrum (PSD)")
plt.show()
Source Power Spectrum (PSD)

Total running time of the script: (0 minutes 6.855 seconds)

Estimated memory usage: 9 MB

Gallery generated by Sphinx-Gallery