Note

Go to the end to download the full example code.

Compute evoked ERS source power using DICS, LCMV beamformer, and dSPM#

Here we examine 3 ways of localizing event-related synchronization (ERS) of beta band activity in this dataset: Somatosensory using DICS, LCMV beamformer, and dSPM applied to active and baseline covariance matrices.

# Authors: Luke Bloy <luke.bloy@gmail.com>
#          Eric Larson <larson.eric.d@gmail.com>
#          Alex Rockhill <aprockhill@mailbox.org>
#
# License: BSD-3-Clause

import mne
import numpy as np
from mne.beamformer import apply_dics_csd, apply_lcmv_cov, make_dics, make_lcmv
from mne.cov import compute_covariance
from mne.datasets import somato
from mne.minimum_norm import apply_inverse_cov, make_inverse_operator
from mne.time_frequency import csd_tfr

import mne_gui_addons

print(__doc__)

Reading the raw data and creating epochs:

data_path = somato.data_path()
subject = "01"
subjects_dir = data_path / "derivatives" / "freesurfer" / "subjects"
task = "somato"
raw_fname = data_path / f"sub-{subject}" / "meg" / f"sub-{subject}_task-{task}_meg.fif"

# crop to 5 minutes to save memory
raw = mne.io.read_raw_fif(raw_fname).crop(0, 300)

# We are interested in the beta band (12-30 Hz)
raw.load_data().filter(12, 30)

# The DICS beamformer currently only supports a single sensor type.
# We'll use the gradiometers in this example.
picks = mne.pick_types(raw.info, meg="grad", exclude="bads")

# Read epochs
events = mne.find_events(raw)
epochs = mne.Epochs(
    raw, events, event_id=1, tmin=-1.5, tmax=2, picks=picks, preload=True, decim=3
)

# Read forward operator and point to freesurfer subject directory
fwd_fname = (
    data_path / "derivatives" / f"sub-{subject}" / f"sub-{subject}_task-{task}-fwd.fif"
)
fwd = mne.read_forward_solution(fwd_fname)

Opening raw data file /home/circleci/mne_data/MNE-somato-data/sub-01/meg/sub-01_task-somato_meg.fif...
    Range : 237600 ... 506999 =    791.189 ...  1688.266 secs
Ready.
Reading 0 ... 90092  =      0.000 ...   299.999 secs...
Filtering raw data in 1 contiguous segment
Setting up band-pass filter from 12 - 30 Hz

FIR filter parameters
---------------------
Designing a one-pass, zero-phase, non-causal bandpass filter:
- Windowed time-domain design (firwin) method
- Hamming window with 0.0194 passband ripple and 53 dB stopband attenuation
- Lower passband edge: 12.00
- Lower transition bandwidth: 3.00 Hz (-6 dB cutoff frequency: 10.50 Hz)
- Upper passband edge: 30.00 Hz
- Upper transition bandwidth: 7.50 Hz (-6 dB cutoff frequency: 33.75 Hz)
- Filter length: 331 samples (1.102 s)

[Parallel(n_jobs=1)]: Done  17 tasks      | elapsed:    0.1s
[Parallel(n_jobs=1)]: Done  71 tasks      | elapsed:    0.3s
[Parallel(n_jobs=1)]: Done 161 tasks      | elapsed:    0.6s
[Parallel(n_jobs=1)]: Done 287 tasks      | elapsed:    1.2s
37 events found on stim channel STI 014
Event IDs: [1]
Not setting metadata
37 matching events found
Setting baseline interval to [-1.498464098687909, 0.0] s
Applying baseline correction (mode: mean)
0 projection items activated
Using data from preloaded Raw for 37 events and 1052 original time points (prior to decimation) ...
0 bad epochs dropped
Reading forward solution from /home/circleci/mne_data/MNE-somato-data/derivatives/sub-01/sub-01_task-somato-fwd.fif...
    Reading a source space...
    [done]
    Reading a source space...
    [done]
    2 source spaces read
    Desired named matrix (kind = 3523 (FIFF_MNE_FORWARD_SOLUTION_GRAD)) not available
    Read MEG forward solution (8155 sources, 306 channels, free orientations)
    Source spaces transformed to the forward solution coordinate frame

Compute covariances and cross-spectral density#

ERS activity starts at 0.5 seconds after stimulus onset. Because these data have been processed by MaxFilter directly (rather than MNE-Python’s version), we have to be careful to compute the rank with a more conservative threshold in order to get the correct data rank (64). Once this is used in combination with an advanced covariance estimator like “shrunk”, the rank will be correctly preserved.

rank = mne.compute_rank(epochs, tol=1e-6, tol_kind="relative")
win_active = (0.5, 1.5)
win_baseline = (-1, 0)
cov_baseline = compute_covariance(
    epochs,
    tmin=win_baseline[0],
    tmax=win_baseline[1],
    method="shrunk",
    rank=rank,
    verbose=True,
)
cov_active = compute_covariance(
    epochs,
    tmin=win_active[0],
    tmax=win_active[1],
    method="shrunk",
    rank=rank,
    verbose=True,
)

# when the covariance objects are added together, they are scaled by the size
# of the window used to create them so that the average is properly weighted
cov_common = cov_baseline + cov_active
cov_baseline.plot(epochs.info)

freqs = np.logspace(np.log10(12), np.log10(30), 9)

# time-frequency decomposition
epochs_tfr = mne.time_frequency.tfr_morlet(
    epochs,
    freqs=freqs,
    n_cycles=freqs / 2,
    return_itc=False,
    average=False,
    output="complex",
)
epochs_tfr.decimate(20)  # decimate for speed

# compute cross-spectral density matrices
csd = csd_tfr(epochs_tfr, tmin=-1, tmax=1.5)
csd_baseline = csd_tfr(epochs_tfr, tmin=win_baseline[0], tmax=win_baseline[1])
csd_ers = csd_tfr(epochs_tfr, tmin=win_active[0], tmax=win_active[1])

csd_baseline.plot()
  • Gradiometers covariance
  • Gradiometers covariance
  • Cross-spectral density, 12.0 Hz., 13.5 Hz., 15.1 Hz., 16.9 Hz., 19.0 Hz., 21.3 Hz., 23.9 Hz., 26.8 Hz., 30.0 Hz.
Computing rank from data with rank=None
    Estimated rank (grad): 64
    GRAD: rank 64 computed from 204 data channels with 0 projectors
    Setting small GRAD eigenvalues to zero (without PCA)
Reducing data rank from 204 -> 64
Estimating covariance using SHRUNK
Done.
Number of samples used : 3737
[done]
    Setting small GRAD eigenvalues to zero (without PCA)
Reducing data rank from 204 -> 64
Estimating covariance using SHRUNK
Done.
Number of samples used : 3737
[done]
Computing rank from covariance with rank=None
    Using tolerance 2.5e-13 (2.2e-16 eps * 204 dim * 5.6  max singular value)
    Estimated rank (grad): 64
    GRAD: rank 64 computed from 204 data channels with 0 projectors
NOTE: tfr_morlet() is a legacy function. New code should use .compute_tfr(method="morlet").
[Parallel(n_jobs=1)]: Done  17 tasks      | elapsed:    0.2s
[Parallel(n_jobs=1)]: Done  71 tasks      | elapsed:    1.0s
[Parallel(n_jobs=1)]: Done 161 tasks      | elapsed:    2.1s

Compute some source estimates#

Here we will use DICS, LCMV beamformer, and dSPM.

See Compute source power using DICS beamformer for more information about DICS.

def _gen_dics(csd, ers_csd, csd_baseline, fwd):
    filters = make_dics(
        epochs.info,
        fwd,
        csd.mean(),
        pick_ori="max-power",
        reduce_rank=True,
        real_filter=True,
        rank=rank,
    )
    stc_base, freqs = apply_dics_csd(csd_baseline.mean(), filters)
    stc_act, freqs = apply_dics_csd(csd_ers.mean(), filters)
    stc_act /= stc_base
    return stc_act


# generate lcmv source estimate
def _gen_lcmv(active_cov, cov_baseline, common_cov, fwd):
    filters = make_lcmv(
        epochs.info, fwd, common_cov, reg=0.05, noise_cov=None, pick_ori="max-power"
    )
    stc_base = apply_lcmv_cov(cov_baseline, filters)
    stc_act = apply_lcmv_cov(cov_active, filters)
    stc_act /= stc_base
    return stc_act


# generate mne/dSPM source estimate
def _gen_mne(cov_active, cov_baseline, cov_common, fwd, info, method="dSPM"):
    inverse_operator = make_inverse_operator(info, fwd, cov_common)
    stc_act = apply_inverse_cov(
        cov_active, info, inverse_operator, method=method, verbose=True
    )
    stc_base = apply_inverse_cov(
        cov_baseline, info, inverse_operator, method=method, verbose=True
    )
    stc_act /= stc_base
    return stc_act


# Compute source estimates
stc_dics = _gen_dics(csd, csd_ers, csd_baseline, fwd)
stc_lcmv = _gen_lcmv(cov_active, cov_baseline, cov_common, fwd)
stc_dspm = _gen_mne(cov_active, cov_baseline, cov_common, fwd, epochs.info)

Identifying common channels ...
Dropped the following channels:
['MEG 0731', 'MEG 1841', 'MEG 1341', 'MEG 0111', 'MEG 2411', 'MEG 2231', 'MEG 1111', 'MEG 0631', 'MEG 2241', 'MEG 2521', 'MEG 1521', 'MEG 1241', 'MEG 1141', 'MEG 2141', 'MEG 0431', 'MEG 1531', 'MEG 0621', 'MEG 1211', 'MEG 2041', 'MEG 2621', 'MEG 1831', 'MEG 2541', 'MEG 0331', 'MEG 1731', 'MEG 1441', 'MEG 0221', 'MEG 0541', 'MEG 1641', 'MEG 2341', 'MEG 2121', 'MEG 2011', 'MEG 0531', 'MEG 0411', 'MEG 2321', 'MEG 1541', 'MEG 1911', 'MEG 1721', 'MEG 1311', 'MEG 2111', 'MEG 2031', 'MEG 0141', 'MEG 1511', 'MEG 2131', 'MEG 0441', 'MEG 1011', 'MEG 1621', 'MEG 0511', 'MEG 1021', 'MEG 2641', 'MEG 2511', 'MEG 1811', 'MEG 0321', 'MEG 2531', 'MEG 1041', 'MEG 0231', 'MEG 1331', 'MEG 0121', 'MEG 1431', 'MEG 1711', 'MEG 0341', 'MEG 0721', 'MEG 0941', 'MEG 2631', 'MEG 1121', 'MEG 2021', 'MEG 1131', 'MEG 2331', 'MEG 0241', 'MEG 2221', 'MEG 0741', 'MEG 1741', 'MEG 0211', 'MEG 0421', 'MEG 0521', 'MEG 0311', 'MEG 0911', 'MEG 2421', 'MEG 2441', 'MEG 1221', 'MEG 1031', 'MEG 1821', 'MEG 0921', 'MEG 1421', 'MEG 0811', 'MEG 2431', 'MEG 2611', 'MEG 1941', 'MEG 1321', 'MEG 1631', 'MEG 0711', 'MEG 1411', 'MEG 0131', 'MEG 2311', 'MEG 1231', 'MEG 2211', 'MEG 0821', 'MEG 1931', 'MEG 1611', 'MEG 1921', 'MEG 0641', 'MEG 0611', 'MEG 0931']
Computing inverse operator with 204 channels.
    204 out of 204 channels remain after picking
Selected 204 channels
Creating the depth weighting matrix...
Whitening the forward solution.
Computing rank from covariance with rank={'grad': 64}
    Setting small GRAD eigenvalues to zero (without PCA)
Creating the source covariance matrix
Adjusting source covariance matrix.
Computing rank from covariance with rank={'grad': 64}
Computing DICS spatial filters...
Computing beamformer filters for 8155 sources
Filter computation complete
Computing DICS source power...
[done]
Computing DICS source power...
[done]
Computing rank from covariance with rank='info'
    GRAD: rank 64 after 0 projectors applied to 204 channels
Computing rank from covariance with rank='info'
    GRAD: rank 64 after 0 projectors applied to 204 channels
Making LCMV beamformer with rank {'grad': 64}
Computing inverse operator with 204 channels.
    204 out of 306 channels remain after picking
Selected 204 channels
Whitening the forward solution.
Computing rank from covariance with rank={'grad': 64}
    Setting small GRAD eigenvalues to zero (without PCA)
Creating the source covariance matrix
Adjusting source covariance matrix.
Computing beamformer filters for 8155 sources
Filter computation complete
Converting forward solution to surface orientation
    No patch info available. The standard source space normals will be employed in the rotation to the local surface coordinates....
    Converting to surface-based source orientations...
    [done]
Computing inverse operator with 204 channels.
    204 out of 306 channels remain after picking
Selected 204 channels
Creating the depth weighting matrix...
    204 planar channels
    limit = 7615/8155 = 10.004172
    scale = 5.17919e-08 exp = 0.8
Applying loose dipole orientations to surface source spaces: 0.2
Whitening the forward solution.
Computing rank from covariance with rank=None
    Using tolerance 2.6e-13 (2.2e-16 eps * 204 dim * 5.8  max singular value)
    Estimated rank (grad): 64
    GRAD: rank 64 computed from 204 data channels with 0 projectors
    Setting small GRAD eigenvalues to zero (without PCA)
Creating the source covariance matrix
Adjusting source covariance matrix.
Computing SVD of whitened and weighted lead field matrix.
    largest singular value = 2.37501
    scaling factor to adjust the trace = 9.82112e+18 (nchan = 204 nzero = 140)
Preparing the inverse operator for use...
    Scaled noise and source covariance from nave = 1 to nave = 1
    Created the regularized inverter
    The projection vectors do not apply to these channels.
    Created the whitener using a noise covariance matrix with rank 64 (140 small eigenvalues omitted)
    Computing noise-normalization factors (dSPM)...
[done]
Applying inverse operator to "cov"...
    Picked 204 channels from the data
    Computing inverse...
    Eigenleads need to be weighted ...
    Computing residual...
    Explained  92.5% variance
    dSPM...
    Combining the current components...
[done]
Preparing the inverse operator for use...
    Scaled noise and source covariance from nave = 1 to nave = 1
    Created the regularized inverter
    The projection vectors do not apply to these channels.
    Created the whitener using a noise covariance matrix with rank 64 (140 small eigenvalues omitted)
    Computing noise-normalization factors (dSPM)...
[done]
Applying inverse operator to "cov"...
    Picked 204 channels from the data
    Computing inverse...
    Eigenleads need to be weighted ...
    Computing residual...
    Explained  92.5% variance
    dSPM...
    Combining the current components...
[done]

Plot source estimates#

DICS:

brain_dics = stc_dics.plot(
    hemi="rh",
    subjects_dir=subjects_dir,
    subject=subject,
    time_label="DICS source power in the 12-30 Hz frequency band",
)
evoked ers source power
Using control points [1.48510357 1.63084332 3.11482832]

LCMV:

brain_lcmv = stc_lcmv.plot(
    hemi="rh",
    subjects_dir=subjects_dir,
    subject=subject,
    time_label="LCMV source power in the 12-30 Hz frequency band",
)
evoked ers source power
Using control points [1.47975676 1.58609544 3.05115173]

dSPM:

brain_dspm = stc_dspm.plot(
    hemi="rh",
    subjects_dir=subjects_dir,
    subject=subject,
    time_label="dSPM source power in the 12-30 Hz frequency band",
)
evoked ers source power
Using control points [1.71153425 1.88773084 3.99261235]

Use volume source estimate with time-frequency resolution#

# make a volume source space
surface = subjects_dir / subject / "bem" / "inner_skull.surf"
vol_src = mne.setup_volume_source_space(
    subject=subject,
    subjects_dir=subjects_dir,
    surface=surface,
    pos=10,
    add_interpolator=False,
)  # just for speed!

conductivity = (0.3,)  # one layer for MEG
model = mne.make_bem_model(
    subject=subject,
    ico=3,  # just for speed
    conductivity=conductivity,
    subjects_dir=subjects_dir,
)
bem = mne.make_bem_solution(model)

trans = fwd["info"]["mri_head_t"]
vol_fwd = mne.make_forward_solution(
    raw.info,
    trans=trans,
    src=vol_src,
    bem=bem,
    meg=True,
    eeg=True,
    mindist=5.0,
    n_jobs=1,
    verbose=True,
)

# Compute source estimate using MNE solver
snr = 3.0
lambda2 = 1.0 / snr**2
method = "MNE"  # use MNE method (could also be dSPM or sLORETA)

# make a different inverse operator for each frequency so as to properly
# whiten the sensor data
inverse_operator = list()
for freq_idx in range(epochs_tfr.freqs.size):
    # for each frequency, compute a separate covariance matrix
    cov_baseline = csd_baseline.get_data(index=freq_idx, as_cov=True)
    cov_baseline["data"] = cov_baseline["data"].real  # only normalize by real
    # then use that covariance matrix as normalization for the inverse
    # operator
    inverse_operator.append(
        mne.minimum_norm.make_inverse_operator(epochs.info, vol_fwd, cov_baseline)
    )

# finally, compute the stcs for each epoch and frequency
stcs = mne.minimum_norm.apply_inverse_tfr_epochs(
    epochs_tfr, inverse_operator, lambda2, method=method, pick_ori="vector"
)

Boundary surface file : /home/circleci/mne_data/MNE-somato-data/derivatives/freesurfer/subjects/01/bem/inner_skull.surf
grid                  : 10.0 mm
mindist               : 5.0 mm
MRI volume            : /home/circleci/mne_data/MNE-somato-data/derivatives/freesurfer/subjects/01/mri/T1.mgz

Reading /home/circleci/mne_data/MNE-somato-data/derivatives/freesurfer/subjects/01/mri/T1.mgz...

Loaded bounding surface from /home/circleci/mne_data/MNE-somato-data/derivatives/freesurfer/subjects/01/bem/inner_skull.surf (10242 nodes)
Surface CM = (   0.1  -21.3   32.9) mm
Surface fits inside a sphere with radius   97.6 mm
Surface extent:
    x =  -76.0 ...   77.4 mm
    y = -105.3 ...   75.3 mm
    z =  -48.1 ...  103.3 mm
Grid extent:
    x =  -80.0 ...   80.0 mm
    y = -110.0 ...   80.0 mm
    z =  -50.0 ...  110.0 mm
5780 sources before omitting any.
3719 sources after omitting infeasible sources not within 0.0 - 97.6 mm.
Source spaces are in MRI coordinates.
Checking that the sources are inside the surface and at least    5.0 mm away (will take a few...)
Checking surface interior status for 3719 points...
    Found  525/3719 points inside  an interior sphere of radius   50.2 mm
    Found    0/3719 points outside an exterior sphere of radius   97.6 mm
    Found 1568/3194 points outside using surface Qhull
    Found  100/1626 points outside using solid angles
    Total 2051/3719 points inside the surface
Interior check completed in 6767.7 ms
    1668 source space points omitted because they are outside the inner skull surface.
    403 source space points omitted because of the    5.0-mm distance limit.
1648 sources remaining after excluding the sources outside the surface and less than    5.0 mm inside.
Adjusting the neighborhood info.
Source space : MRI voxel -> MRI (surface RAS)
    0.010000 0.000000 0.000000     -80.00 mm
    0.000000 0.010000 0.000000    -110.00 mm
    0.000000 0.000000 0.010000     -50.00 mm
    0.000000 0.000000 0.000000       1.00
MRI volume : MRI voxel -> MRI (surface RAS)
    -0.001000 0.000000 0.000000     128.00 mm
    0.000000 0.000000 0.001000    -128.00 mm
    0.000000 -0.001000 0.000000     128.00 mm
    0.000000 0.000000 0.000000       1.00
MRI volume : MRI (surface RAS) -> RAS (non-zero origin)
    1.000000 0.000000 0.000000       2.93 mm
    0.000000 1.000000 0.000000      49.85 mm
    0.000000 0.000000 1.000000     -24.58 mm
    0.000000 0.000000 0.000000       1.00
Creating the BEM geometry...
Going from 5th to 3th subdivision of an icosahedron (n_tri: 20480 -> 1280)
inner skull CM is   0.15 -21.35  32.94 mm
Surfaces passed the basic topology checks.
Complete.

Homogeneous model surface loaded.
Computing the linear collocation solution...
    Matrix coefficients...
        inner skull (642) -> inner skull (642) ...
    Inverting the coefficient matrix...
Solution ready.
BEM geometry computations complete.
Source space          : <SourceSpaces: [<volume, shape=(np.int64(17), np.int64(20), np.int64(17)), n_used=1648>] MRI (surface RAS) coords, subject '01', ~1.4 MB>
MRI -> head transform : instance of Transform
Measurement data      : instance of Info
Conductor model   : instance of ConductorModel
Accurate field computations
Do computations in head coordinates
Free source orientations

Read 1 source spaces a total of 1648 active source locations

Coordinate transformation: MRI (surface RAS) -> head
    0.998865 0.023957 0.041176       0.49 mm
    -0.035225 0.953344 0.299821      24.87 mm
    -0.032072 -0.300931 0.953106      17.96 mm
    0.000000 0.000000 0.000000       1.00

Read 306 MEG channels from info
105 coil definitions read
Coordinate transformation: MEG device -> head
    0.997640 -0.061033 -0.031446       3.42 mm
    0.063766 0.993449 0.094829       6.64 mm
    0.025452 -0.096610 0.994997      54.42 mm
    0.000000 0.000000 0.000000       1.00
MEG coil definitions created in head coordinates.
Source spaces are now in head coordinates.

Employing the head->MRI coordinate transform with the BEM model.
BEM model instance of ConductorModel is now set up

Source spaces are in head coordinates.
Checking that the sources are inside the surface and at least    5.0 mm away (will take a few...)
Checking surface interior status for 1648 points...
    Found  531/1648 points inside  an interior sphere of radius   50.7 mm
    Found    0/1648 points outside an exterior sphere of radius   97.3 mm
    Found    0/1117 points outside using surface Qhull
    Found    0/1117 points outside using solid angles
    Total 1648/1648 points inside the surface
Interior check completed in 310.5 ms

Checking surface interior status for 306 points...
    Found   0/306 points inside  an interior sphere of radius   50.7 mm
    Found 306/306 points outside an exterior sphere of radius   97.3 mm
    Found   0/  0 points outside using surface Qhull
    Found   0/  0 points outside using solid angles
    Total 0/306 points inside the surface
Interior check completed in 50.0 ms

Composing the field computation matrix...
Computing MEG at 1648 source locations (free orientations)...

Finished.
Computing inverse operator with 204 channels.
    204 out of 306 channels remain after picking
Selected 204 channels
Creating the depth weighting matrix...
    204 planar channels
    limit = 1559/1648 = 10.054558
    scale = 4.43241e-08 exp = 0.8
Whitening the forward solution.
Computing rank from covariance with rank=None
    Using tolerance 5.1e-13 (2.2e-16 eps * 204 dim * 11  max singular value)
    Estimated rank (grad): 64
    GRAD: rank 64 computed from 204 data channels with 0 projectors
    Setting small GRAD eigenvalues to zero (without PCA)
Creating the source covariance matrix
Adjusting source covariance matrix.
Computing SVD of whitened and weighted lead field matrix.
    largest singular value = 2.44896
    scaling factor to adjust the trace = 4.03799e+18 (nchan = 204 nzero = 140)
Computing inverse operator with 204 channels.
    204 out of 306 channels remain after picking
Selected 204 channels
Creating the depth weighting matrix...
    204 planar channels
    limit = 1559/1648 = 10.054558
    scale = 4.43241e-08 exp = 0.8
Whitening the forward solution.
Computing rank from covariance with rank=None
    Using tolerance 5.8e-13 (2.2e-16 eps * 204 dim * 13  max singular value)
    Estimated rank (grad): 64
    GRAD: rank 64 computed from 204 data channels with 0 projectors
    Setting small GRAD eigenvalues to zero (without PCA)
Creating the source covariance matrix
Adjusting source covariance matrix.
Computing SVD of whitened and weighted lead field matrix.
    largest singular value = 2.44948
    scaling factor to adjust the trace = 3.38374e+18 (nchan = 204 nzero = 140)
Computing inverse operator with 204 channels.
    204 out of 306 channels remain after picking
Selected 204 channels
Creating the depth weighting matrix...
    204 planar channels
    limit = 1559/1648 = 10.054558
    scale = 4.43241e-08 exp = 0.8
Whitening the forward solution.
Computing rank from covariance with rank=None
    Using tolerance 5.7e-13 (2.2e-16 eps * 204 dim * 13  max singular value)
    Estimated rank (grad): 64
    GRAD: rank 64 computed from 204 data channels with 0 projectors
    Setting small GRAD eigenvalues to zero (without PCA)
Creating the source covariance matrix
Adjusting source covariance matrix.
Computing SVD of whitened and weighted lead field matrix.
    largest singular value = 2.39781
    scaling factor to adjust the trace = 3.57564e+18 (nchan = 204 nzero = 140)
Computing inverse operator with 204 channels.
    204 out of 306 channels remain after picking
Selected 204 channels
Creating the depth weighting matrix...
    204 planar channels
    limit = 1559/1648 = 10.054558
    scale = 4.43241e-08 exp = 0.8
Whitening the forward solution.
Computing rank from covariance with rank=None
    Using tolerance 4.5e-13 (2.2e-16 eps * 204 dim * 9.9  max singular value)
    Estimated rank (grad): 64
    GRAD: rank 64 computed from 204 data channels with 0 projectors
    Setting small GRAD eigenvalues to zero (without PCA)
Creating the source covariance matrix
Adjusting source covariance matrix.
Computing SVD of whitened and weighted lead field matrix.
    largest singular value = 2.45123
    scaling factor to adjust the trace = 3.97057e+18 (nchan = 204 nzero = 140)
Computing inverse operator with 204 channels.
    204 out of 306 channels remain after picking
Selected 204 channels
Creating the depth weighting matrix...
    204 planar channels
    limit = 1559/1648 = 10.054558
    scale = 4.43241e-08 exp = 0.8
Whitening the forward solution.
Computing rank from covariance with rank=None
    Using tolerance 3e-13 (2.2e-16 eps * 204 dim * 6.7  max singular value)
    Estimated rank (grad): 64
    GRAD: rank 64 computed from 204 data channels with 0 projectors
    Setting small GRAD eigenvalues to zero (without PCA)
Creating the source covariance matrix
Adjusting source covariance matrix.
Computing SVD of whitened and weighted lead field matrix.
    largest singular value = 2.49264
    scaling factor to adjust the trace = 4.81685e+18 (nchan = 204 nzero = 140)
Computing inverse operator with 204 channels.
    204 out of 306 channels remain after picking
Selected 204 channels
Creating the depth weighting matrix...
    204 planar channels
    limit = 1559/1648 = 10.054558
    scale = 4.43241e-08 exp = 0.8
Whitening the forward solution.
Computing rank from covariance with rank=None
    Using tolerance 1.6e-13 (2.2e-16 eps * 204 dim * 3.5  max singular value)
    Estimated rank (grad): 64
    GRAD: rank 64 computed from 204 data channels with 0 projectors
    Setting small GRAD eigenvalues to zero (without PCA)
Creating the source covariance matrix
Adjusting source covariance matrix.
Computing SVD of whitened and weighted lead field matrix.
    largest singular value = 2.43023
    scaling factor to adjust the trace = 6.16806e+18 (nchan = 204 nzero = 140)
Computing inverse operator with 204 channels.
    204 out of 306 channels remain after picking
Selected 204 channels
Creating the depth weighting matrix...
    204 planar channels
    limit = 1559/1648 = 10.054558
    scale = 4.43241e-08 exp = 0.8
Whitening the forward solution.
Computing rank from covariance with rank=None
    Using tolerance 9.7e-14 (2.2e-16 eps * 204 dim * 2.2  max singular value)
    Estimated rank (grad): 64
    GRAD: rank 64 computed from 204 data channels with 0 projectors
    Setting small GRAD eigenvalues to zero (without PCA)
Creating the source covariance matrix
Adjusting source covariance matrix.
Computing SVD of whitened and weighted lead field matrix.
    largest singular value = 2.23101
    scaling factor to adjust the trace = 7.9978e+18 (nchan = 204 nzero = 140)
Computing inverse operator with 204 channels.
    204 out of 306 channels remain after picking
Selected 204 channels
Creating the depth weighting matrix...
    204 planar channels
    limit = 1559/1648 = 10.054558
    scale = 4.43241e-08 exp = 0.8
Whitening the forward solution.
Computing rank from covariance with rank=None
    Using tolerance 6.8e-14 (2.2e-16 eps * 204 dim * 1.5  max singular value)
    Estimated rank (grad): 64
    GRAD: rank 64 computed from 204 data channels with 0 projectors
    Setting small GRAD eigenvalues to zero (without PCA)
Creating the source covariance matrix
Adjusting source covariance matrix.
Computing SVD of whitened and weighted lead field matrix.
    largest singular value = 2.26033
    scaling factor to adjust the trace = 1.07631e+19 (nchan = 204 nzero = 140)
Computing inverse operator with 204 channels.
    204 out of 306 channels remain after picking
Selected 204 channels
Creating the depth weighting matrix...
    204 planar channels
    limit = 1559/1648 = 10.054558
    scale = 4.43241e-08 exp = 0.8
Whitening the forward solution.
Computing rank from covariance with rank=None
    Using tolerance 4.6e-14 (2.2e-16 eps * 204 dim * 1  max singular value)
    Estimated rank (grad): 64
    GRAD: rank 64 computed from 204 data channels with 0 projectors
    Setting small GRAD eigenvalues to zero (without PCA)
Creating the source covariance matrix
Adjusting source covariance matrix.
Computing SVD of whitened and weighted lead field matrix.
    largest singular value = 2.28605
    scaling factor to adjust the trace = 1.62809e+19 (nchan = 204 nzero = 140)
Not setting metadata
37 matching events found
No baseline correction applied
0 projection items activated
Preparing the inverse operator for use...
    Scaled noise and source covariance from nave = 1 to nave = 1
    Created the regularized inverter
    The projection vectors do not apply to these channels.
    Created the whitener using a noise covariance matrix with rank 64 (140 small eigenvalues omitted)
Picked 204 channels from the data
Computing inverse...
    Eigenleads need to be weighted ...
Processing epoch : 1 / 37
Processing epoch : 2 / 37
Processing epoch : 3 / 37
Processing epoch : 4 / 37
Processing epoch : 5 / 37
Processing epoch : 6 / 37
Processing epoch : 7 / 37
Processing epoch : 8 / 37
Processing epoch : 9 / 37
Processing epoch : 10 / 37
Processing epoch : 11 / 37
Processing epoch : 12 / 37
Processing epoch : 13 / 37
Processing epoch : 14 / 37
Processing epoch : 15 / 37
Processing epoch : 16 / 37
Processing epoch : 17 / 37
Processing epoch : 18 / 37
Processing epoch : 19 / 37
Processing epoch : 20 / 37
Processing epoch : 21 / 37
Processing epoch : 22 / 37
Processing epoch : 23 / 37
Processing epoch : 24 / 37
Processing epoch : 25 / 37
Processing epoch : 26 / 37
Processing epoch : 27 / 37
Processing epoch : 28 / 37
Processing epoch : 29 / 37
Processing epoch : 30 / 37
Processing epoch : 31 / 37
Processing epoch : 32 / 37
Processing epoch : 33 / 37
Processing epoch : 34 / 37
Processing epoch : 35 / 37
Processing epoch : 36 / 37
Processing epoch : 37 / 37
[done]
Not setting metadata
37 matching events found
No baseline correction applied
0 projection items activated
Preparing the inverse operator for use...
    Scaled noise and source covariance from nave = 1 to nave = 1
    Created the regularized inverter
    The projection vectors do not apply to these channels.
    Created the whitener using a noise covariance matrix with rank 64 (140 small eigenvalues omitted)
Picked 204 channels from the data
Computing inverse...
    Eigenleads need to be weighted ...
Processing epoch : 1 / 37
Processing epoch : 2 / 37
Processing epoch : 3 / 37
Processing epoch : 4 / 37
Processing epoch : 5 / 37
Processing epoch : 6 / 37
Processing epoch : 7 / 37
Processing epoch : 8 / 37
Processing epoch : 9 / 37
Processing epoch : 10 / 37
Processing epoch : 11 / 37
Processing epoch : 12 / 37
Processing epoch : 13 / 37
Processing epoch : 14 / 37
Processing epoch : 15 / 37
Processing epoch : 16 / 37
Processing epoch : 17 / 37
Processing epoch : 18 / 37
Processing epoch : 19 / 37
Processing epoch : 20 / 37
Processing epoch : 21 / 37
Processing epoch : 22 / 37
Processing epoch : 23 / 37
Processing epoch : 24 / 37
Processing epoch : 25 / 37
Processing epoch : 26 / 37
Processing epoch : 27 / 37
Processing epoch : 28 / 37
Processing epoch : 29 / 37
Processing epoch : 30 / 37
Processing epoch : 31 / 37
Processing epoch : 32 / 37
Processing epoch : 33 / 37
Processing epoch : 34 / 37
Processing epoch : 35 / 37
Processing epoch : 36 / 37
Processing epoch : 37 / 37
[done]
Not setting metadata
37 matching events found
No baseline correction applied
0 projection items activated
Preparing the inverse operator for use...
    Scaled noise and source covariance from nave = 1 to nave = 1
    Created the regularized inverter
    The projection vectors do not apply to these channels.
    Created the whitener using a noise covariance matrix with rank 64 (140 small eigenvalues omitted)
Picked 204 channels from the data
Computing inverse...
    Eigenleads need to be weighted ...
Processing epoch : 1 / 37
Processing epoch : 2 / 37
Processing epoch : 3 / 37
Processing epoch : 4 / 37
Processing epoch : 5 / 37
Processing epoch : 6 / 37
Processing epoch : 7 / 37
Processing epoch : 8 / 37
Processing epoch : 9 / 37
Processing epoch : 10 / 37
Processing epoch : 11 / 37
Processing epoch : 12 / 37
Processing epoch : 13 / 37
Processing epoch : 14 / 37
Processing epoch : 15 / 37
Processing epoch : 16 / 37
Processing epoch : 17 / 37
Processing epoch : 18 / 37
Processing epoch : 19 / 37
Processing epoch : 20 / 37
Processing epoch : 21 / 37
Processing epoch : 22 / 37
Processing epoch : 23 / 37
Processing epoch : 24 / 37
Processing epoch : 25 / 37
Processing epoch : 26 / 37
Processing epoch : 27 / 37
Processing epoch : 28 / 37
Processing epoch : 29 / 37
Processing epoch : 30 / 37
Processing epoch : 31 / 37
Processing epoch : 32 / 37
Processing epoch : 33 / 37
Processing epoch : 34 / 37
Processing epoch : 35 / 37
Processing epoch : 36 / 37
Processing epoch : 37 / 37
[done]
Not setting metadata
37 matching events found
No baseline correction applied
0 projection items activated
Preparing the inverse operator for use...
    Scaled noise and source covariance from nave = 1 to nave = 1
    Created the regularized inverter
    The projection vectors do not apply to these channels.
    Created the whitener using a noise covariance matrix with rank 64 (140 small eigenvalues omitted)
Picked 204 channels from the data
Computing inverse...
    Eigenleads need to be weighted ...
Processing epoch : 1 / 37
Processing epoch : 2 / 37
Processing epoch : 3 / 37
Processing epoch : 4 / 37
Processing epoch : 5 / 37
Processing epoch : 6 / 37
Processing epoch : 7 / 37
Processing epoch : 8 / 37
Processing epoch : 9 / 37
Processing epoch : 10 / 37
Processing epoch : 11 / 37
Processing epoch : 12 / 37
Processing epoch : 13 / 37
Processing epoch : 14 / 37
Processing epoch : 15 / 37
Processing epoch : 16 / 37
Processing epoch : 17 / 37
Processing epoch : 18 / 37
Processing epoch : 19 / 37
Processing epoch : 20 / 37
Processing epoch : 21 / 37
Processing epoch : 22 / 37
Processing epoch : 23 / 37
Processing epoch : 24 / 37
Processing epoch : 25 / 37
Processing epoch : 26 / 37
Processing epoch : 27 / 37
Processing epoch : 28 / 37
Processing epoch : 29 / 37
Processing epoch : 30 / 37
Processing epoch : 31 / 37
Processing epoch : 32 / 37
Processing epoch : 33 / 37
Processing epoch : 34 / 37
Processing epoch : 35 / 37
Processing epoch : 36 / 37
Processing epoch : 37 / 37
[done]
Not setting metadata
37 matching events found
No baseline correction applied
0 projection items activated
Preparing the inverse operator for use...
    Scaled noise and source covariance from nave = 1 to nave = 1
    Created the regularized inverter
    The projection vectors do not apply to these channels.
    Created the whitener using a noise covariance matrix with rank 64 (140 small eigenvalues omitted)
Picked 204 channels from the data
Computing inverse...
    Eigenleads need to be weighted ...
Processing epoch : 1 / 37
Processing epoch : 2 / 37
Processing epoch : 3 / 37
Processing epoch : 4 / 37
Processing epoch : 5 / 37
Processing epoch : 6 / 37
Processing epoch : 7 / 37
Processing epoch : 8 / 37
Processing epoch : 9 / 37
Processing epoch : 10 / 37
Processing epoch : 11 / 37
Processing epoch : 12 / 37
Processing epoch : 13 / 37
Processing epoch : 14 / 37
Processing epoch : 15 / 37
Processing epoch : 16 / 37
Processing epoch : 17 / 37
Processing epoch : 18 / 37
Processing epoch : 19 / 37
Processing epoch : 20 / 37
Processing epoch : 21 / 37
Processing epoch : 22 / 37
Processing epoch : 23 / 37
Processing epoch : 24 / 37
Processing epoch : 25 / 37
Processing epoch : 26 / 37
Processing epoch : 27 / 37
Processing epoch : 28 / 37
Processing epoch : 29 / 37
Processing epoch : 30 / 37
Processing epoch : 31 / 37
Processing epoch : 32 / 37
Processing epoch : 33 / 37
Processing epoch : 34 / 37
Processing epoch : 35 / 37
Processing epoch : 36 / 37
Processing epoch : 37 / 37
[done]
Not setting metadata
37 matching events found
No baseline correction applied
0 projection items activated
Preparing the inverse operator for use...
    Scaled noise and source covariance from nave = 1 to nave = 1
    Created the regularized inverter
    The projection vectors do not apply to these channels.
    Created the whitener using a noise covariance matrix with rank 64 (140 small eigenvalues omitted)
Picked 204 channels from the data
Computing inverse...
    Eigenleads need to be weighted ...
Processing epoch : 1 / 37
Processing epoch : 2 / 37
Processing epoch : 3 / 37
Processing epoch : 4 / 37
Processing epoch : 5 / 37
Processing epoch : 6 / 37
Processing epoch : 7 / 37
Processing epoch : 8 / 37
Processing epoch : 9 / 37
Processing epoch : 10 / 37
Processing epoch : 11 / 37
Processing epoch : 12 / 37
Processing epoch : 13 / 37
Processing epoch : 14 / 37
Processing epoch : 15 / 37
Processing epoch : 16 / 37
Processing epoch : 17 / 37
Processing epoch : 18 / 37
Processing epoch : 19 / 37
Processing epoch : 20 / 37
Processing epoch : 21 / 37
Processing epoch : 22 / 37
Processing epoch : 23 / 37
Processing epoch : 24 / 37
Processing epoch : 25 / 37
Processing epoch : 26 / 37
Processing epoch : 27 / 37
Processing epoch : 28 / 37
Processing epoch : 29 / 37
Processing epoch : 30 / 37
Processing epoch : 31 / 37
Processing epoch : 32 / 37
Processing epoch : 33 / 37
Processing epoch : 34 / 37
Processing epoch : 35 / 37
Processing epoch : 36 / 37
Processing epoch : 37 / 37
[done]
Not setting metadata
37 matching events found
No baseline correction applied
0 projection items activated
Preparing the inverse operator for use...
    Scaled noise and source covariance from nave = 1 to nave = 1
    Created the regularized inverter
    The projection vectors do not apply to these channels.
    Created the whitener using a noise covariance matrix with rank 64 (140 small eigenvalues omitted)
Picked 204 channels from the data
Computing inverse...
    Eigenleads need to be weighted ...
Processing epoch : 1 / 37
Processing epoch : 2 / 37
Processing epoch : 3 / 37
Processing epoch : 4 / 37
Processing epoch : 5 / 37
Processing epoch : 6 / 37
Processing epoch : 7 / 37
Processing epoch : 8 / 37
Processing epoch : 9 / 37
Processing epoch : 10 / 37
Processing epoch : 11 / 37
Processing epoch : 12 / 37
Processing epoch : 13 / 37
Processing epoch : 14 / 37
Processing epoch : 15 / 37
Processing epoch : 16 / 37
Processing epoch : 17 / 37
Processing epoch : 18 / 37
Processing epoch : 19 / 37
Processing epoch : 20 / 37
Processing epoch : 21 / 37
Processing epoch : 22 / 37
Processing epoch : 23 / 37
Processing epoch : 24 / 37
Processing epoch : 25 / 37
Processing epoch : 26 / 37
Processing epoch : 27 / 37
Processing epoch : 28 / 37
Processing epoch : 29 / 37
Processing epoch : 30 / 37
Processing epoch : 31 / 37
Processing epoch : 32 / 37
Processing epoch : 33 / 37
Processing epoch : 34 / 37
Processing epoch : 35 / 37
Processing epoch : 36 / 37
Processing epoch : 37 / 37
[done]
Not setting metadata
37 matching events found
No baseline correction applied
0 projection items activated
Preparing the inverse operator for use...
    Scaled noise and source covariance from nave = 1 to nave = 1
    Created the regularized inverter
    The projection vectors do not apply to these channels.
    Created the whitener using a noise covariance matrix with rank 64 (140 small eigenvalues omitted)
Picked 204 channels from the data
Computing inverse...
    Eigenleads need to be weighted ...
Processing epoch : 1 / 37
Processing epoch : 2 / 37
Processing epoch : 3 / 37
Processing epoch : 4 / 37
Processing epoch : 5 / 37
Processing epoch : 6 / 37
Processing epoch : 7 / 37
Processing epoch : 8 / 37
Processing epoch : 9 / 37
Processing epoch : 10 / 37
Processing epoch : 11 / 37
Processing epoch : 12 / 37
Processing epoch : 13 / 37
Processing epoch : 14 / 37
Processing epoch : 15 / 37
Processing epoch : 16 / 37
Processing epoch : 17 / 37
Processing epoch : 18 / 37
Processing epoch : 19 / 37
Processing epoch : 20 / 37
Processing epoch : 21 / 37
Processing epoch : 22 / 37
Processing epoch : 23 / 37
Processing epoch : 24 / 37
Processing epoch : 25 / 37
Processing epoch : 26 / 37
Processing epoch : 27 / 37
Processing epoch : 28 / 37
Processing epoch : 29 / 37
Processing epoch : 30 / 37
Processing epoch : 31 / 37
Processing epoch : 32 / 37
Processing epoch : 33 / 37
Processing epoch : 34 / 37
Processing epoch : 35 / 37
Processing epoch : 36 / 37
Processing epoch : 37 / 37
[done]
Not setting metadata
37 matching events found
No baseline correction applied
0 projection items activated
Preparing the inverse operator for use...
    Scaled noise and source covariance from nave = 1 to nave = 1
    Created the regularized inverter
    The projection vectors do not apply to these channels.
    Created the whitener using a noise covariance matrix with rank 64 (140 small eigenvalues omitted)
Picked 204 channels from the data
Computing inverse...
    Eigenleads need to be weighted ...
Processing epoch : 1 / 37
Processing epoch : 2 / 37
Processing epoch : 3 / 37
Processing epoch : 4 / 37
Processing epoch : 5 / 37
Processing epoch : 6 / 37
Processing epoch : 7 / 37
Processing epoch : 8 / 37
Processing epoch : 9 / 37
Processing epoch : 10 / 37
Processing epoch : 11 / 37
Processing epoch : 12 / 37
Processing epoch : 13 / 37
Processing epoch : 14 / 37
Processing epoch : 15 / 37
Processing epoch : 16 / 37
Processing epoch : 17 / 37
Processing epoch : 18 / 37
Processing epoch : 19 / 37
Processing epoch : 20 / 37
Processing epoch : 21 / 37
Processing epoch : 22 / 37
Processing epoch : 23 / 37
Processing epoch : 24 / 37
Processing epoch : 25 / 37
Processing epoch : 26 / 37
Processing epoch : 27 / 37
Processing epoch : 28 / 37
Processing epoch : 29 / 37
Processing epoch : 30 / 37
Processing epoch : 31 / 37
Processing epoch : 32 / 37
Processing epoch : 33 / 37
Processing epoch : 34 / 37
Processing epoch : 35 / 37
Processing epoch : 36 / 37
Processing epoch : 37 / 37
[done]

Plot volume source estimates#

viewer = mne_gui_addons.view_vol_stc(
    stcs, subject=subject, subjects_dir=subjects_dir, src=vol_src, inst=epochs_tfr
)
viewer.go_to_extreme()  # show the maximum intensity source vertex
viewer.set_cmap(vmin=0.25, vmid=0.8)
viewer.set_3d_view(azimuth=40, elevation=35, distance=350)
GUI

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

Estimated memory usage: 3287 MB

Gallery generated by Sphinx-Gallery